Directional Subsurface Scattering Jeppe Revall Frisvad June 2020 - - PowerPoint PPT Presentation

Directional Subsurface Scattering

Jeppe Revall Frisvad June 2020

Frisvad, J. R., Hachisuka, T., and Kjeldsen, T. K. Directional dipole model for subsurface scattering. ACM Transactions on Graphics 34(1),

• pp. 5:1-5:12, November 2014. Presented at SIGGRAPH 2015.

Materials (scattering and absorption of light)

◮ Optical properties (index of refraction, n(λ) = n′(λ) + i n′′(λ)). ◮ Reflectance distribution functions, S(xi, ωi; xo, ωo).

xi xo n1 n2

BSSRDF ◮ The BSSRDF (Bidirectional Scattering-Surface Reflectance Distribution Function) describes surface and subsurface scattering.

An estimator for subsurface scattering

◮ The reflected radiance equation: Lr(xo, ωo) =

• A
• 2π

S(xi, ωi; xo, ωo)Li(xi, ωi) cos θi dωi dAi . ◮ Monte Carlo estimator: Lr,N,M(xo, ωo) = 1 NM

M

• p=1

N

• q=1

S(xi,p, ωi,q; xo, ωo)Li(xi,p, ωi,q) cos θi pdf(xi,p)pdf( ωi,q) . ◮ Common direction sampling pdf (cosine-weighted hemisphere): pdf( ωi,q) = ωi,q · ni π = cos θi π . ◮ Common area sampling pdf (triangle mesh): pdf(xi,p) = pdf(△)pdf(xi,p,△) = A△ Aℓ 1 A△ = 1 Aℓ .

Sampling a cosine-weighted hemisphere (ambient occlusion)

◮ Material: S(xi, ωi; xo, ωo) = ρd(xo) π δ(xo − xi) . ◮ Sampler: pdf( ωi,q) = ωi,q · ni π = cos θi π . ◮ Estimator: Lr,N(xo, ωo) = 1 N

N

• q=1

ρd(xo) π Li(xo, ωi,q) cos θi pdf( ωi,q) = ρd(xo) 1 N

M

• q=1

Li(xo, ωi,q) .

Sampling a triangle mesh (area lights, soft shadows)

◮ Material: S(xi, ωi; xo, ωo) = ρd(xo) π δ(xo − xi) . ◮ Sampler:

• ωi,q =

xℓ,q − xi xℓ,q − xi pdf(xℓ,q) = pdf(△)pdf(xℓ,q,△) = A△ Aℓ 1 A△ . ◮ Estimator: Lr,N(xo, ωo) = ρd(xo) π 1 N

N

• q=1

Le(xℓ,q, − ωi,q)V (xℓ,q, xo) cos θi cos θℓ xℓ,q − xi2 Aℓ .

Sampling for subsurface scattering

◮ Material: S(xi, ωi; xo, ωo) = . . . . ◮ Sampler: pdf( ωi,q) = ωi,q · ni π = cos θi π . pdf(xi,p) = pdf(△)pdf(xi,p,△) = A△ Aℓ 1 A△ = 1 Aℓ .

• ωo
• ωi,2

xo

• no

xi,2

• ni,2

Ai,2 xi,1 dωi,1

• ni,1

Ai,1

References

• Frisvad, J. R., Hachisuka, T., and Kjeldsen, T. K. Directional dipole model for subsurface scattering. ACM Transactions on Graphics 34(1),
• pp. 5:1–5:12, November 2014. Presented at SIGGRAPH 2015.
• Dal Corso, A., and Frisvad, J. R. Point cloud method for rendering BSSRDFs. Technical Report, Technical University of Denmark, 2018.

Splitting up the BSSRDF

◮ Bidirectional Scattering-Surface Reflectance Distribution Function: S = S(xi, ωi; xo, ωo) . ◮ Away from sources and boundaries, we can use diffusion. ◮ Splitting up the BSSRDF S = T12(S(0) + S(1) + Sd)T21 . where

◮ T12 and T21 are Fresnel transmittance terms (using ωi, ωo). ◮ S(0) is the direct transmission part (using Dirac δ-functions). ◮ S(1) is the single scattering part (using all arguments). ◮ Sd is the diffusive part (multiple scattering, using |xo − xi|).

◮ We distribute the single scattering to the other terms using the delta-Eddington approximation: S = T12(SδE + Sd)T21 , and generalize the model such that Sd = Sd(xi, ωi; xo).

Analytical models for subsurface scattering

standard dipole S(xi, ωi; xo, ωo) = T12( ωi)(S1 + Sd(xo − xi))T21( ωo) . directional dipole S(xi, ωi; xo, ωo) = T12( ωi)(SδE + Sd(xi, ωi; xo))T21( ωo) . ◮ Directions ( ωi, ωo) also require surface normals ( ni, no) to get angles (θi, θo). ◮ T12 and T21 are Fresnel transmittances. ◮ S1 and SδE are fully directional (depend on xi, ωi, xo, ωo, and normals).

Diffusion theory

◮ Think of multiple scattering as a diffusion process. ◮ In diffusion theory, we use quantities that describe the light field in an element of volume of the scattering medium. ◮ Total flux, or fluence, is defined by φ(x) =

• 4π

L(x, ω) dω .

x y z x dx dy dz

◮ We find an expression for φ by solving the diffusion equation (D∇2 − σa)φ(x) = −q(x) + 3D ∇·Q(x) , where σa and D are absorption and diffusion coefficients, while q and Q are zeroth and first order source terms.

Deriving a BSSRDF

◮ Assume that emerging light is diffuse due to a large number of scattering events: Sd(xi, ωi; xo, ωo) = Sd(xi, ωi; xo). ◮ Integrating emerging diffuse radiance over outgoing directions, we find Sd = Cφ(η) φ − CE(η) D no·∇φ Φ 4πCφ(1/η) , where

◮ Φ is the flux entering the medium at xi. ◮ no is the surface normal at the point of emergence xo. ◮ Cφ and CE depend on the relative index of refraction η and are polynomial fits of different hemispherical integrals of the Fresnel transmittance.

◮ This connects the BSSRDF and the diffusion theory. ◮ To get an analytical model, we use a special case solution for the diffusion equation (an expression for φ). ◮ Then, “all” we need to do is to find ∇φ (do the math) and deal with boundary conditions (build a plausible model).

Point source diffusion or ray source diffusion

standard dipole ◮ Point source diffusion [Bothe 1941; 1942] φ(r) =

Φ 4πD e−σtrr r

, where r = |xo − xi| and σtr =

• σa/D is the effective

transport coefficient. directional dipole ◮ Ray source diffusion [Menon et al. 2005a; 2005b] φ(r, θ) =

Φ 4πD e−σtrr r

• 1 + 3D 1+σtrr

r

cos θ

• ,

where θ is the angle between the refracted ray and xo − xi.

Diffusive part of the standard dipole

Sd(r) = α′ 4π2 zr(1 + σtrdr)e−σtrdr d3

r

+ zv(1 + σtrdv)e−σtrdv d3

v

• .

◮ Distances:

◮ zr = Λ . ◮ zv = Λ + 4AD . ◮ dr(r) =

• r 2 + z2

r .

◮ dv(r) =

• r 2 + z2

v .

◮ Optical properties (η = n2/n1, σs, σa, g):

◮ Reduced scattering coefficient: σ′

s = σs(1 − g) .

◮ Reduced extinction coefficient: σ′

t = σ′ s + σa .

◮ Reduced scattering albedo: α′ = σ′

s/σ′ t .

◮ Transport mean free path: Λ = 1/σ′

t .

◮ Diffusion coefficient: D = Λ/3 . ◮ Transport coefficient: σtr =

• σa/D .

◮ Reflection parameter: A(η) (ratio of polynomial fits).

r z-axis z

• zv

r ωi ωo n dv z = 0 dr real source virtual source

• bserver

incident light n2 n1

Directional subsurface scattering when disregarding the boundary

S′

d(x,

ω12, r) =

1 4Cφ(1/η) 1 4π2 e−σtrr r 3

• Cφ(η)
• r 2

D + 3(1 + σtrr) x ·

ω12

• − CE(η)
• 3D(1 + σtrr)

ω12 · no −

• (1 + σtrr) + 3D 3(1+σtrr)+(σtrr)2

r 2

x · ω12

• x ·

no

• ,

where Cφ(η) and CE(η) are polynomial fits. ◮ Additional dependencies:

◮ Normal: no . ◮ Optical properites: η, D, σtr .

◮ Note the exponential term: e−σtrr . ◮ Normal incidence: ω12 · no = ±1 . ◮ Plane (half-space): x · no ≈ 0 . ◮ normal incidence on plane: x · ω12 ≈ 0 . ◮ r → xo − xi for xo − xi → ∞ .

d

Dipole configuration (method of mirror images)

d

◮ We place the “real” ray source at the boundary and reflect it in an extrapolated boundary to place the “virtual” ray source. ◮ Distance to the extrapolated boundary [Davison 1958]: de = 2.131 D/

• 1 − 3Dσa .

◮ In case of a refractive boundary (η1 = η2), the distance is Ade with A = 1 − CE(η) 2Cφ(η) .

Modified tangent plane

d

◮ The dipole assumes a semi-infinite medium. ◮ We assume that the boundary contains the vector xo − xi and that it is perpendicular to the plane spanned by ni and xo − xi. ◮ The normal of the assumed boundary plane is then

• n ∗

i = xo − xi

|xo − xi| × ni × (xo − xi) | ni × (xo − xi)|,

• r

n ∗

i =

ni if xo = xi. and the virtual source is given by xv = xi + 2Ade n ∗

i , dv = |xv − xi| ,

ωv = ω12 − 2( ω12 · n ∗

i )

n ∗

i .

Distance to the real source (handling the singularity)

r z-axis z

• zv

r ωi ωo n dv z = 0 dr real source virtual source

• bserver

incident light n2 n1

standard dipole dr =

• r2 + z2

r .

d

directional dipole dr = r ? ◮ Emergent radiance is an integral over z of a Hankel transform of a Green function which is Fourier transformed in x and y. ◮ Approximate analytic evaluation is possible if r is corrected to R2 = r2 + (z′ + de)2 . ◮ The resulting model for z′ = 0 corresponds to the standard dipole where z′ = zr and de is replaced by the virtual source.

Distance to the real source (handling the singularity)

• r2 + d2

e

de xi xo

• ω12

θ θ0

• no

− no r β

◮ Since we neither have normal incidence nor xo in the tangent plane, we modify the distance correction: R2 = r2 + z′2 + d2

e − 2z′de cos β .

◮ The integral over z can be reformulated as an integral along the refracted ray. ◮ We can approximate this integral by choosing an offset D∗ along the refracted ray. Then z′ = D∗| cos θ0|.

Diffusive part of the directional dipole

Sd(xi, ωi; xo) = S′

d(xo − xi,

ω12, dr) − S′

d(xo − xv,

ωv, dv) , ◮ Real source:

◮ ω12 = η−1(( ωi · ni) ni − ωi) − ni

• 1 − η−2(1 − (

ωi · ni)2) . ◮ d2

r =

xo − xi2 + Dµ0(Dµ0 − 2de cos β) for µ0 > 0 xo − xi2 + 1/(3σt)2

• therwise ,

with µ0 = cos θ0 = − no · ω12 and cos β = −

• r 2−(x·ω12)2

r 2+d2

e

.

◮ Virtual source:

◮ Modified normal:

• n ∗

i = xo−xi xo−xi ×

• ni×(xo−xi)
• ni×(xo−xi),
• r

n ∗

i =

ni if xo = xi. ◮ xv = xi + 2Ade n ∗

i , dv = |xv − xi| ,

ωv = ω12 − 2( ω12 · n ∗

i )

n ∗

i .

d

Rejection control

◮ The exponential attenuation e−σtrd appears in all analytical BSSRDFs and d → xo − xi for xo − xi → ∞ . ◮ We should exploit this. ◮ Russian roulette:

sample ξ ∈ [0, 1] uniformly; if (ξ < P1)

call event 1; divide by p1;

else if (ξ < P2)

call event 2; divide by p2;

else if (ξ < P3)

. . .

else if (ξ < P4)

. . .

◮ When sampling xi, use Russian roulette with p1(xi) = P1(xi) = e−σtrxo−xi to accept or reject a sample.

1 1 2 3 4

Discrete pdf

p - probability 1 1 2 3 4

Discrete cdf

P - cumulave

Progressive rendering of subsurface scattering

◮ The equation for reflected radiance: Lr(xo, ωo) =

• A
• 2π

S(xi, ωi; xo, ωo)Li(xi, ωi) cos θi dωi dAi . ◮ Initialize a frame by storing samples of transmitted light. ◮ For each sample:

◮ Sample a point (xi,p) on the surface of the scattering material by sampling a random triangle and then a random point in the triangle: pdf(xi,p) = 1/Aℓ . ◮ Sample a ray direction ωi,q using a cosine-weighted hemisphere and trace they ray to collect incident light Li: pdf( ωi,q) = cos θi/π . ◮ Use ωi,q to find the direction of the transmitted/refracted ray and the Fresnel transmittance T12. ◮ Store the transmitted radiance: Lt = T12Li cos θi pdf(xi,p)pdf( ωi,q) = T12LiπAℓ.

Progressive rendering of subsurface scattering

◮ The equation for reflected radiance: Lr(xo, ωo) =

• A
• 2π

S(xi, ωi; xo, ωo)Li(xi, ωi) cos θi dωi dAi . ◮ For a ray hitting xo with direction − ωo:

◮ Compute Fresnel transmittance T21 of the ray refracting from inside to ωo. ◮ Loop through the NM samples using exponential distance attenuation as the probability of acceptance in a Russian roulette (rejection control). ◮ Use Lt of accepted samples and T21 together with the analytical expression for S to Monte Carlo integrate the rendering equation. ◮ The Monte Carlo estimator for the diffusive part is (we use N = 1): Ld,N,M(xo, ωo) = 1 NM

M

• p=1

N

• q=1

S(xi,p, ωi,q; xo, ωo)Li(xi,p, ωi,q) cos θi pdf(xi,p)pdf( ωi,q) = T21 NM

M

• p=1

N

• q=1

Sd(xi,p, ωi,q; xo, ωo)T12Li(xi,p, ωi,q)πAℓ e−σtrxo−xi,p

• ξ < e−σtrxo−xi,p

.

Profiles (diffuse reflectance curves)

R (x)

d

• 16
• 12
• 8
• 4

4 8 12 16 10 1 10 0 10 -1 10 -2 10 -3 10 -4 10 -5

dipole btpole qntzd

• urs

ptrace

dipole btpole qntzd

• urs

ptrace

◮ The directional model comes closer than the diffuse analytical models to measured and simulated diffuse reflectance curves.