Skip to content

fix(bsdl): Correct colored Fresnel sampling #2158

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
tdavidovicNV wants to merge 1 commit into
base: main
Choose a base branch
from
+142 −23
Open

fix(bsdl): Correct colored Fresnel sampling #2158

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
11 changes: 10 additions & 1 deletion src/libbsdl/CMakeLists.txt
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -7,4 +7,13 @@ include(bsdl.cmake)
# BSDL has small tables for doing spectral render in sRGB
# but for other color spaces tell this function which ones
# and it will bake Jakob-Hanika coefficient tables.
add_bsdl_library(BSDL) # SPECTRAL_COLOR_SPACES "ACEScg")
add_bsdl_library(BSDL) # SPECTRAL_COLOR_SPACES "ACEScg")


if (OSL_BUILD_TESTS AND BUILD_TESTING)
add_executable (bsdl_test bsdl_test.cpp)
set_target_properties (bsdl_test PROPERTIES FOLDER "Unit Tests")
target_include_directories (bsdl_test BEFORE PRIVATE ${OpenImageIO_INCLUDES})
target_link_libraries (bsdl_test PRIVATE BSDL OpenImageIO::OpenImageIO)
add_test (unit_bsdl ${CMAKE_RUNTIME_OUTPUT_DIRECTORY}/bsdl_test)
endif()
93 changes: 93 additions & 0 deletions src/libbsdl/bsdl_test.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,93 @@
// Copyright Contributors to the Open Shading Language project.
// SPDX-License-Identifier: BSD-3-Clause
// https://github.com/AcademySoftwareFoundation/OpenShadingLanguage

#define BSDL_UNROLL()

#include <BSDL/config.h>

using BSDLConfig = bsdl::BSDLDefaultConfig;

#include <BSDL/MTX/bsdf_schlick_impl.h>

#include <OpenImageIO/unittest.h>

using namespace bsdl;



struct TestDielectricBSDF : mtx::DielectricBSDF<mtx::SchlickFresnel> {
using Base = mtx::DielectricBSDF<mtx::SchlickFresnel>;
using Base::Base;
using Base::reflection_probability;
};



static void
test_colored_schlick_sampling()
{
const float rgb[] = { 0.1f, 0.2f, 0.8f };
const Power F([&](int i) { return rgb[i]; }, 0.0f);
const mtx::SchlickFresnel fresnel(F, F, 5.0f, 1.5f, false);
const GGXDist dist(0.25f, 0.0f);
const Imath::V3f wo(0.0f, 0.0f, 1.0f);
const TestDielectricBSDF bsdf(dist, fresnel, wo.z, 0.25f, true, 0.0f);

const float probability = (rgb[0] + rgb[1] + rgb[2]) / 3.0f;
OIIO_CHECK_EQUAL_THRESH(bsdf.reflection_probability(F), probability, 1e-6f);
OIIO_CHECK_EQUAL_THRESH(bsdf.reflection_probability(Power::UNIT()), 1.0f,
1e-6f);

const Sample reflected = bsdf.sample(wo, 0.5f, 0.5f, probability - 0.01f);
const Sample transmitted = bsdf.sample(wo, 0.5f, 0.5f, probability + 0.01f);
OIIO_CHECK_ASSERT(reflected.wi.z > 0.0f);
OIIO_CHECK_ASSERT(transmitted.wi.z < 0.0f);

const float rgb2[] = { 0.2f, 0.4f, 0.6f };
const Power F2([&](int i) { return rgb2[i]; }, 0.0f);
const mtx::SchlickFresnel fresnel2(F2, F2, 5.0f, 1.5f, false);
const TestDielectricBSDF bsdf2(dist, fresnel2, wo.z, 0.25f, true, 0.0f);
const float probability2 = (rgb2[0] + rgb2[1] + rgb2[2]) / 3.0f;

const Sample reflected2 = bsdf2.sample(wo, 0.5f, 0.5f, 0.3f);
const Sample transmitted2 = bsdf2.sample(wo, 0.5f, 0.5f, 0.8f);
const Sample reflected1 = bsdf.sample(wo, 0.5f, 0.5f, 0.3f);
const Sample transmitted1 = bsdf.sample(wo, 0.5f, 0.5f, 0.8f);
OIIO_CHECK_EQUAL_THRESH(reflected1.pdf / reflected2.pdf,
probability / probability2, 1e-6f);
const float transmission_probability = 1.0f - probability;
const float transmission_probability2 = 1.0f - probability2;
OIIO_CHECK_EQUAL_THRESH(transmitted1.pdf / transmitted2.pdf,
transmission_probability
/ transmission_probability2,
1e-6f);
for (int i = 0; i < 3; ++i) {
OIIO_CHECK_EQUAL_THRESH(reflected1.weight[i] / reflected2.weight[i],
(rgb[i] * probability2)
/ (rgb2[i] * probability),
1e-6f);
OIIO_CHECK_EQUAL_THRESH(transmitted1.weight[i] / transmitted2.weight[i],
((1.0f - rgb[i]) * transmission_probability2)
/ ((1.0f - rgb2[i])
* transmission_probability),
1e-6f);
}

const float spectral[] = { 0.2f, 0.4f, 0.6f, 0.8f };
const Power Fs([&](int i) { return spectral[i]; }, 500.0f);
const mtx::SchlickFresnel spectral_fresnel(Fs, Fs, 5.0f, 1.5f, false);
const TestDielectricBSDF spectral_bsdf(dist, spectral_fresnel, wo.z, 0.25f,
true, 500.0f);
OIIO_CHECK_EQUAL_THRESH(spectral_bsdf.reflection_probability(Fs), 0.5f,
1e-6f);
}



int
main(int /*argc*/, char* /*argv*/[])
{
test_colored_schlick_sampling();
return unit_test_failures;
}
5 changes: 4 additions & 1 deletion src/libbsdl/include/BSDL/MTX/bsdf_dielectric_decl.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -42,7 +42,7 @@ template<typename Fresnel> struct DielectricBSDF {

BSDL_INLINE_METHOD
DielectricBSDF(const GGXDist& dist, const Fresnel& fresnel, float cosNO,
float roughness, bool dorefr);
float roughness, bool dorefr, float lambda_0);

DielectricBSDF() = default;

Expand All @@ -58,9 +58,12 @@ template<typename Fresnel> struct DielectricBSDF {
static constexpr const char* NS = "mtx";

protected:
BSDL_INLINE_METHOD float reflection_probability(const Power& F) const;

GGXDist d;
Fresnel f;
float E_ms;
float lambda_0;
bool dorefr;
};

Expand Down
53 changes: 33 additions & 20 deletions src/libbsdl/include/BSDL/MTX/bsdf_dielectric_impl.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -67,8 +67,9 @@ template<typename Fresnel>
BSDL_INLINE_METHOD
DielectricBSDF<Fresnel>::DielectricBSDF(const GGXDist& dist,
const Fresnel& fresnel, float cosNO,
float roughness, bool dorefr)
: d(dist), f(fresnel), dorefr(dorefr)
float roughness, bool dorefr,
float lambda_0)
: d(dist), f(fresnel), lambda_0(lambda_0), dorefr(dorefr)
{
if (!dorefr) {
TabulatedEnergyCurve<spi::MiniMicrofacetGGX> curve(roughness, 0.0f);
Expand All @@ -83,7 +84,7 @@ DielectricReflFront::DielectricReflFront(float cosNO, float roughness_index,
: DielectricBSDF<DielectricFresnel>(
GGXDist(roughness_index, 0),
DielectricFresnel::from_table_index(fresnel_index, false), cosNO,
roughness_index, false)
roughness_index, false, 1)
{
}

Expand All @@ -93,7 +94,7 @@ DielectricBothFront::DielectricBothFront(float cosNO, float roughness_index,
: DielectricBSDF<DielectricFresnel>(
GGXDist(roughness_index, 0),
DielectricFresnel::from_table_index(fresnel_index, false), cosNO,
roughness_index, true)
roughness_index, true, 1)

@aconty aconty Sep 3, 2026

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Don't we need to pass lambda_0 here too?

@tdavidovicNV tdavidovicNV Sep 3, 2026

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

As far as I understand, these two only exist for LUT generation, use scalar DielectricFresnel where all 4 lanes have the same value and therefore the lamba_0 doesn't really matter. But I agree that the 1 isn't the most self-explanatory here.

We could push this whole thing into the Fresnel calculation, so the BSDF wouldn't need to know about the lambda, only the colored SchlickFresnel would, but that would move the responsibility of calculating the probabilities from the BSDF (where it lives now) to Fresnel, which would get a reflection_probability interface. I tried to make the change as minimal as possible, but it is definitely a tradeoff to consider.

@aconty aconty Sep 3, 2026

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I remember now. As a sanity check, can you check if the generated luts changed? Would it be better to use 0?

@tdavidovicNV tdavidovicNV Sep 4, 2026

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@aconty so with lambda_0 == 1, the results are bit identical. I have just tried with = 0 and it is actually worse.

Before the proposed change, we just used F.max(), which picked on of the identical lanes. With lambda_0 = 1, we sum all 4 lanes and then divide by 4, and get the identical number. With lambda_0 = 0, we would sum up only first 3 and divide by 3 and, by the magic of powers of 2 being easier to optimize, we get a different result.

Given the standard IEEE754 approach of having Guard, Round, and Sticky bits in the FPUs, you pretty much have 2 bits of extra precision in the arithmetic, so it makes perfect sense that exactly 4x bigger/smaller makes no difference (you don't even hit the Sticky bit).

{
}

Expand All @@ -103,10 +104,19 @@ DielectricBothBack::DielectricBothBack(float cosNO, float roughness_index,
: DielectricBSDF<DielectricFresnel>(
GGXDist(roughness_index, 0),
DielectricFresnel::from_table_index(fresnel_index, true), cosNO,
roughness_index, true)
roughness_index, true, 1)

@aconty aconty Sep 3, 2026

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

and here?

{
}

template<typename Fresnel>
BSDL_INLINE_METHOD float
DielectricBSDF<Fresnel>::reflection_probability(const Power& F) const
{
// Fresnel may return Power::UNIT(), so mask the inactive RGB lane before
// averaging.
return F.cliped_rgb(lambda_0).avg(lambda_0);
}

template<typename Fresnel>
BSDL_INLINE_METHOD Sample
DielectricBSDF<Fresnel>::eval(Imath::V3f wo, Imath::V3f wi) const
Expand All @@ -125,18 +135,19 @@ DielectricBSDF<Fresnel>::eval(Imath::V3f wo, Imath::V3f wi) const
const float D = d.D(m);
const float G1 = d.G1(wo);
const Power F = f.eval(cosMO);
if (F.max() <= 0)
const float P = reflection_probability(F);
if (P <= 0)
return {};
if constexpr (BSDLConfig::use_bvn_refraction) {
// Reflection optimized density
const float D_refl_D = d.D_refl_D(wo, m);
const float D_refl = D_refl_D * D;
const Power out = F * (d.G2_G1(wi, wo) * G1 / (D_refl_D * F.max()));
const float pdf = D_refl / (4.0f * cosNO) * F.max();
const Power out = F * (d.G2_G1(wi, wo) * G1 / (D_refl_D * P));
const float pdf = D_refl / (4.0f * cosNO) * P;
return { wi, out, pdf, 0 };
} else {
const Power out = F * d.G2_G1(wi, wo);
const float pdf = (G1 * D * F.max()) / (4.0f * cosNO);
const Power out = F * (d.G2_G1(wi, wo) / P);
const float pdf = (G1 * D * P) / (4.0f * cosNO);
return { wi, out, pdf, 0 };
}
} else if (cosNI < 0) {
Expand All @@ -148,8 +159,10 @@ DielectricBSDF<Fresnel>::eval(Imath::V3f wo, Imath::V3f wi) const
const float cosHI = Ht.dot(wi);
if (cosHO <= 0 || cosHI >= 0)
return {};
const Power Ft = Power::UNIT() - f.eval(cosHO);
if (Ht.z <= 0 || Ft.max() <= 0)
const Power F = f.eval(cosHO);
const Power Ft = Power::UNIT() - F;
const float Pt = 1 - reflection_probability(F);
if (Ht.z <= 0 || Pt <= 0)
return {};
const float D = d.D(Ht);
const float G1 = d.G1(wo);
Expand All @@ -159,16 +172,15 @@ DielectricBSDF<Fresnel>::eval(Imath::V3f wo, Imath::V3f wi) const
// Reflection optimized density
const float D_refl_D = d.D_refl_D(wo, Ht);
const float D_refl = D_refl_D * D;
float pdf = D_refl * J * Ft.max();
float pdf = D_refl * J * Pt;
const Power out = Ft
* (d.G2_G1({ wi.x, wi.y, -wi.z }, wo) * G1
/ (D_refl_D * Ft.max()));
/ (D_refl_D * Pt));
return { wi, out, pdf, 0 };
} else {
const Power out = Ft
* (d.G2_G1({ wi.x, wi.y, -wi.z }, wo) / Ft.max());
const Power out = Ft * (d.G2_G1({ wi.x, wi.y, -wi.z }, wo) / Pt);

float pdf = J * G1 * D * Ft.max();
float pdf = J * G1 * D * Pt;
return { wi, out, pdf, 0 };
}

Expand Down Expand Up @@ -196,8 +208,8 @@ DielectricBSDF<Fresnel>::sample(Imath::V3f wo, float randu, float randv,
const float cosMO = wo.dot(m);
if (cosMO <= 0)
return {};
const float F = f.eval(cosMO).max();
bool choose_reflect = randw < F;
const float P = reflection_probability(f.eval(cosMO));
bool choose_reflect = randw < P;
const Imath::V3f wi = choose_reflect ? reflect(wo, m)
: refract(wo, m, f.refraction_eta());
if ((choose_reflect && wi.z <= 0) || (!choose_reflect && wi.z >= 0))
Expand Down Expand Up @@ -245,7 +257,8 @@ DielectricLobe<BSDF_ROOT>::DielectricLobe(T* lobe, const BsdfGlobals& globals,
DielectricFresnel fresnel(globals.relative_eta(IOR), globals.backfacing);
E_ms = 0;
spec = DielectricBSDF<DielectricFresnel>(GGXDist(roughness, aniso, rx < ry),
fresnel, cosNO, roughness, dorefr);
fresnel, cosNO, roughness, dorefr,
globals.lambda_0);
if (dorefl && !dorefr) {
E_ms = TabulatedEnergyCurve<DielectricReflFront>(roughness,
fresnel.table_index())
Expand Down
3 changes: 2 additions & 1 deletion src/libbsdl/include/BSDL/MTX/bsdf_schlick_impl.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -73,7 +73,8 @@ SchlickLobe<BSDF_ROOT>::SchlickLobe(T* lobe, const BsdfGlobals& globals,

E_ms = 0;
spec = DielectricBSDF<SchlickFresnel>(GGXDist(roughness, aniso, rx < ry),
fresnel, cosNO, roughness, dorefr);
fresnel, cosNO, roughness, dorefr,
globals.lambda_0);
if (dorefl && !dorefr) {
// Energy compensation reuses the dielectric Fresnel albedo tables,
// which assumes the Schlick curve matches the true dielectric Fresnel.
Expand Down
Loading