22// SPDX-License-Identifier: BSD-3-Clause
33// https://github.com/AcademySoftwareFoundation/OpenShadingLanguage
44
5-
65#pragma once
76
87#include < BSDL/config.h>
@@ -70,6 +69,50 @@ GGXDist::sample(const Imath::V3f& wo, float randu, float randv) const
7069 return Imath::V3f (ax * N.x , ay * N.y , std::max (N.z , 0 .0f )).normalized ();
7170}
7271
72+ BSDL_INLINE_METHOD Imath::V3f
73+ GGXDist::sample_reflection (const Imath::V3f& wo, float randu, float randv) const
74+ {
75+ // From "Bounded VNDF Sampling for Smith–GGX Reflections" Eto - Tokuyoshi.
76+ // This code is taken from listing 1 almost verbatim.
77+ Imath::V3f i_std = Imath::V3f (wo.x * ax, wo.y * ay, wo.z ).normalized ();
78+ // Sample a spherical cap
79+ const float phi = 2 .0f * PI * randu;
80+ const float a = CLAMP (std::min (ax, ay), 0 .0f , 1 .0f ); // Eq. 6
81+ const float s = 1 + sqrtf (SQR (wo.x ) + SQR (wo.y )); // Omit sgn for a <=1
82+ const float a2 = SQR (a);
83+ const float s2 = SQR (s);
84+ const float k = (1 - a2) * s2 / (s2 + a2 * SQR (wo.z )); // Eq. 5
85+ const float b = k * i_std.z ;
86+ const float z = (1 - randv) * (1 + b) - b;
87+ const float sinTheta = sqrtf (CLAMP (1 - SQR (z), 0 .0f , 1 .0f ));
88+ constexpr auto cos = BSDLConfig::Fast::cosf;
89+ constexpr auto sin = BSDLConfig::Fast::sinf;
90+ Imath::V3f o_std = { sinTheta * cos (phi), sinTheta * sin (phi), z };
91+ // Compute the microfacet normal m
92+ Imath::V3f m_std = i_std + o_std;
93+ Imath::V3f m = Imath::V3f (m_std.x * ax, m_std.y * ay, m_std.z ).normalized ();
94+ // Return the reflection vector o
95+ return 2 * wo.dot (m) * m - wo;
96+ }
97+
98+ BSDL_INLINE_METHOD float
99+ GGXDist::pdf_reflection (const Imath::V3f& wo, const Imath::V3f& wi) const
100+ {
101+ // From "Bounded VNDF Sampling for Smith–GGX Reflections" Eto - Tokuyoshi.
102+ // This code is taken from listing 2 almost verbatim.
103+ const Imath::V3f m = (wo + wi).normalized ();
104+ const float ndf = D (m);
105+ const Imath::V2f ai = { wo.x * ax, wo.y * ay };
106+ const float len2 = ai.dot (ai);
107+ const float t = sqrtf (len2 + SQR (wo.z ));
108+ const float a = CLAMP (std::min (ax, ay), 0 .0f , 1 .0f ); // Eq. 6
109+ const float s = 1 + sqrtf (SQR (wo.x ) + SQR (wo.y )); // Omit sgn for a <=1
110+ const float a2 = SQR (a);
111+ const float s2 = SQR (s);
112+ const float k = (1 - a2) * s2 / (s2 + a2 * SQR (wo.z )); // Eq. 5
113+ return ndf / (2 * (k * wo.z + t)); // Eq. 8 * || dm/do ||
114+ }
115+
73116template <typename BSDF>
74117BSDL_INLINE_METHOD float
75118TabulatedEnergyCurve<BSDF>::interpolate_emiss(int i) const
@@ -299,9 +342,9 @@ MicrofacetMS<Fresnel>::computeFmiss() const
299342}
300343
301344// Turquin style microfacet with multiple scattering
302- template <typename Dist, typename Fresnel>
345+ template <typename Fresnel>
303346BSDL_INLINE Sample
304- eval_turquin_microms_reflection (const Dist & dist, const Fresnel& fresnel,
347+ eval_turquin_microms_reflection (const GGXDist & dist, const Fresnel& fresnel,
305348 float E_ms, const Imath::V3f& wo,
306349 const Imath::V3f& wi)
307350{
@@ -315,8 +358,8 @@ eval_turquin_microms_reflection(const Dist& dist, const Fresnel& fresnel,
315358 float cosMO = m.dot (wo);
316359 const float D = dist.D (m);
317360 const float G1 = dist.G1 (wo);
318- const float out = dist.G2_G1 (wi, wo );
319- float s_pdf = (G1 * D) / (4 .0f * cosNO);
361+ float s_pdf = dist.pdf_reflection (wo, wi );
362+ const float out = dist. G2_G1 (wi, wo) * (G1 * D) / (4 .0f * cosNO * s_pdf );
320363 // fresnel term between outgoing direction and microfacet
321364 const Power F = fresnel.eval (cosMO);
322365 // From "Practical multiple scattering compensation for microfacet models" - Emmanuel Turquin
@@ -329,4 +372,85 @@ eval_turquin_microms_reflection(const Dist& dist, const Fresnel& fresnel,
329372 return { wi, O, s_pdf, -1 /* roughness set by caller */ };
330373}
331374
375+ template <typename Fresnel>
376+ BSDL_INLINE Sample
377+ sample_turquin_microms_reflection (const GGXDist& dist, const Fresnel& fresnel,
378+ float E_ms, const Imath::V3f& wo,
379+ const Imath::V3f& rnd)
380+ {
381+ const float cosNO = wo.z ;
382+ if (cosNO <= 0 )
383+ return {};
384+
385+ Imath::V3f wi = dist.sample_reflection (wo, rnd.x , rnd.y );
386+ if (wi.z <= 0 )
387+ return {};
388+
389+ // evaluate brdf on outgoing direction
390+ return eval_turquin_microms_reflection (dist, fresnel, E_ms, wo, wi);
391+ }
392+
393+ // SPI style microfacet with multiple scattering, with a diffuse lobe
394+ template <typename Dist, typename Fresnel>
395+ BSDL_INLINE Sample
396+ eval_spi_microms_reflection (const Dist& dist, const Fresnel& fresnel,
397+ float E_ms, float E_ms_avg, const Imath::V3f& wo,
398+ const Imath::V3f& wi)
399+ {
400+ const float cosNO = wo.z ;
401+ const float cosNI = wi.z ;
402+ if (cosNI <= 0 || cosNO <= 0 )
403+ return {};
404+
405+ // get half vector
406+ Imath::V3f m = (wo + wi).normalized ();
407+ float cosMO = m.dot (wo);
408+ const float D = dist.D (m);
409+ const float G1 = dist.G1 (wo);
410+ const float out = dist.G2_G1 (wi, wo);
411+ float s_pdf = (G1 * D) / (4 .0f * cosNO);
412+ // fresnel term between outgoing direction and microfacet
413+ const Power F = fresnel.eval (cosMO);
414+ const Power O = out * F;
415+
416+ const Power one = Power (1 , 1 );
417+ // Like Turquin we are using F_avg = F instead of an average fresnel, the
418+ // render integrator has an average effect for high roughness.
419+ // const Power F0 = fresnel.F0();
420+ const Power F_avg = F;
421+ // Evaluate multiple scattering lobe, roughness set by caller. The
422+ // 1 / (1 - F_avg * E_ms_avg) term comes from the series summation.
423+ const Power O_ms = SQR (F_avg) * (1 - E_ms_avg) / (one - F_avg * E_ms_avg);
424+
425+ Sample ec = { wi, O_ms, E_ms * cosNI * ONEOVERPI, -1 };
426+ ec.update (O, s_pdf, 1 - E_ms);
427+
428+ return ec;
429+ }
430+
431+ template <typename Dist, typename Fresnel>
432+ BSDL_INLINE Sample
433+ sample_spi_microms_reflection (const Dist& dist, const Fresnel& fresnel,
434+ float E_ms, float E_ms_avg, const Imath::V3f& wo,
435+ const Imath::V3f& rnd)
436+ {
437+ const float cosNO = wo.z ;
438+ if (cosNO <= 0 )
439+ return {};
440+
441+ Imath::V3f wi;
442+ if (rnd.z < E_ms) {
443+ // sample diffuse energy compensation lobe
444+ wi = sample_cos_hemisphere (rnd.x , rnd.y );
445+ } else {
446+ // sample microfacet (half vector)
447+ // generate outgoing direction
448+ wi = reflect (wo, dist.sample (wo, rnd.x , rnd.y ));
449+ if (wi.z <= 0 )
450+ return {};
451+ }
452+ // evaluate brdf on outgoing direction
453+ return eval_spi_microms_reflection (dist, fresnel, E_ms, E_ms_avg, wo, wi);
454+ }
455+
332456BSDL_LEAVE_NAMESPACE
0 commit comments