Thursday, January 22, 2015

BSSRDF Importance Sampling 2 - Dipole Diffusion Approximation101

Super quick review on the dipole diffusion approximation model.
Compare to the classic BRDF integration:
$$L\left ( x_{o},\vec{\omega _{o}} \right )=L_{e}\left ( x_{o},\vec{\omega _{o}} \right )+\int_{\Omega }L_{i}\left ( x_{i},\vec{\omega _{i}} \right )f_{r}\left ( \vec{\omega _{i}},\vec{\omega _{o}} \right )cos\left ( \theta \right )d\omega$$
This is BSSRDF integration:
$$L\left ( x_{o},\vec{\omega _{o}} \right )=\int_{A}\int_{\Omega  }L_{i}\left ( x_{i},\vec{\omega _{i}} \right )S\left (  x_{i},\vec{\omega _{i}}, x_{o},\vec{\omega _{o}} \right )cos\left ( \theta \right )d\omega dA $$

In short words, to evaluate the radiance come out from surface point \(x_{o}\) with direction \(\omega_{o}\), we need to integrate over the whole surface area (the \(\int_{A }dA\) part ), at each surface point \(x_{i}\), we integrate over the hemisphere(the \(\int_{\Omega }d\omega\) part) that for radiance \(L( x_{i},\vec{\omega _{i}}) \), how much portion arrive at point \(x_{o}\) with direction \(\omega_{o}\)(the \(S( x_{i},\vec{\omega _{i}}, x_{o},\vec{\omega _{o}})\) part)

the BSSRDF \(S(  x_{i},\vec{\omega _{i}}, x_{o},\vec{\omega _{o}})\) is a 4 dimensional function, but most of the time it got simplified to radially symmetric and become a 1D distance funtion \( R_{d}(d)\) where d is distance between \( x_{o}\) and \( x_{i}\).

Dipole Diffusion Approximation made assumptions that:
-light distribution  in highly scattering media tends to become isotropic. (My understanding for this part is that when light enter into medium with really high Albedo, after multiple scattering event happens, the light scatter distribution will become more and more uniform, and approach \(1/4\pi\) at last)
-the surface is semi-infinite flat surface since for finite media there is no general way to solve diffusion equation in analytic way.

By putting two point light source above and under the medium(that's why it's called dipole method), we can get the analytic solution for fluence on the semi-infinite flat surface (OK...."we" not including "me" at this moment......I confess...... Donner 06 thesis Chapter 5 contains a detail derivation on this, while I just grab the result equation directly...)

$$R_{d}\left (r\right )=\frac{{\alpha}'}{4\pi}\biggl [z_{r}( \sigma_{tr}d_{r}+1)\frac{e^{-\sigma_{tr}d_{r}}}{ {d_{r}}^{3} } + z_{v}( \sigma_{tr}d_{v}+1)\frac{e^{-\sigma_{tr}d_{v}}}{ {d_{v}}^{3} }\biggr ]$$

After a topological sort on those scary symbols:
\(\sigma_{a}\): absorb coefficient, describe how much radiance got absorbed in medium
\(\sigma_{s}\): scatter coefficient, describe how much radiance got scattered in medium
\(g\): anisotropy parameter (range from -1 to 1, -1 for fully backward scattering, 1 for fully forward scattering, 0 for isotropy)
\({\sigma_{s}}' =\sigma_{s}(1-g) \): reduced scattering coefficient
\({\sigma_{t}}' ={\sigma_{s}}' + \sigma_{a} \): reduced extinction coefficient
\({\alpha}' = {\sigma_{s}}'/{\sigma_{t}}'\): the reduced albedo
\(\sigma_{tr} = \sqrt{3\sigma_{a}{\sigma_{t}}'}\): the effective transport coefficient
\(F_{r}\): the Fresnel formula for the reflection at a dielectric surface
\(F_{dr} = \int_{\Omega }F_{r}(\eta, \vec{n}\cdot \vec{\omega})(\vec{n}\cdot \vec{\omega})d\omega\): the average diffuse Fresnel reflectance
\(A=(1+F_{dr})/(1-F_{dr})\): the coefficient takes into account the effect of internal diffuse reflection
\(D=\frac{1}{3{\sigma_{t}}'}\): the diffusion constant
\(z_{r} = 1/{\sigma_{t}}'\): the distance beneath surface where we put the positive dipole light
\(z_{v} = z_{r} + 4AD\): the distance above surface where we put the negative dipole light
\(d_{r} =\parallel x-x_{r}\parallel\): the distance from x to the positive dipole light
\(d_{v} =\parallel x-x_{v}\parallel\): the distance from x to the negative dipole light

Hopefully this dipole image illustrates the thing a little bit clearer...

At last, we put into account the Fresnel reflection at the boundary for both the incoming light and the outgoing light, we get:
$$S_{d}\left (x_{i}, \vec{\omega_{i}}, x_{o}, \vec{\omega_{o}}\right)=\frac{1}{\pi}F_{r}(\eta, \vec{\omega_{i}})R_{d}(\parallel x_{i}, x_{o} \parallel)F_{r}(\eta, \vec{\omega_{o}})$$

Image taken from GPU Gems3, showing the BSSRDF falloff pretty quick when the distance between \(x_{o}\) and \(x_{i}\) increases, we can also see that red channel tend to resist more before fully extinct, this sorta explains the idea human skin have that red texture underneath

For Fresnel diffuse reflection \(F_{dr}\) , it's an integration that can be approximated with the following code snippet. It was taken from Donner 06 too:

YUUUP.....we got the diffusion approximation part BSSRDF \(S_{d}(x_{i}, \vec{\omega_{i}}, x_{o}, \vec{\omega_{o}})\) now, and it looks pretty mathmatically daunting (again, for me....). So, I encapsulate all these stuff into the BSSRDF black box, and hopefully I don't need to see them most of the time :P The constructor takes \(\sigma_{a}\), \({\sigma_{s}}'\), refraction index \(\eta\), anisotropy parameter \(g\) , the rest of the stuff it can cook indoor. Something looks like this:

OK! That pretty much sums up the TODO list item1, now moves on to item2, get a separate integration for BSSRDF. This is just open another cut in point for some future implementation, pretty straightforward for small side project scale code base:

Done!

One particular thought pop in my mind while reading Jensen01 was: "For the sake of god!!! To just get something radially symmetrically falloff you really need to be that hardcore? I mean....there should be something so much easier symmetrically falloff model you can use right? like Gaussian distribution maybe?" The fun fact is yes! There are already some research about using multiple weight sum of Gaussian filters to approximate BSSRDF and it can be really close to the dipole model. This GPU gems 3 article mentioned the idea and give you the weight sum parameters directly for skin BSSRDF. If you are interested in the full story under the hood, d'Eon11 "A Quantized-Diffusion Model for Rendering Translucent Materials" is a cool ride on it.

Image taken from GPU Gems3 , we can see the 4 Gaussian weight sum is really close to the Dipole model already


2 comments:

  1. Hi Wei-Feng,

    I am a Visual Basic programmer, I don't know much C++.

    I have written a graphics application in VB. It does all the basic stuff and global illumination with photon mapping and final gathering.

    I am trying to implement subsurface scattering now. I have red the Jensen and Donner papers, but I am stuck at how to implement Jensen's diffusion DIPOLE in code.

    I have coded the Rd(r) expression, but I am at a loss on how to integrate the area around the point I am rendering.

    I looked at your code in C++ to get some ideas, but I still need help. I am rendering a point, how to go about integrating the area around it to get sss with Jensen's dipole? Christensen suggests going inside the surface and spherically shoot sampling rays and look at the illumination at the ray intersections. So I have a bunch of points with illumination close to the point I am rendering, how to go from there and apply Jensen's diffuse dipole? Can you help me? I hope you will have time for that, I can be reached at my email radiomilo@sodetel.net.lb. Hope to hear from you! Thank you. Camil.

    ReplyDelete
  2. Hi, to integrate samples around the area. There are two main stream methods: Monte Carlo or density estimation. If you have bunch of points with illumination close to the point you are rendering, you can probably try density estimation. The rough idea is just divide the sum of illumination points (each multiply by Rd(r)) by pi*r^2 to get the result. The dipole integrator in pbrt-v2 is a good example of this kind of density estimation method (you probably need an spatial accelerator like kd-tree for fast neighbor query during render time though):
    https://github.com/mmp/pbrt-v2/blob/master/src/integrators/dipolesubsurface.cpp

    ReplyDelete