# What is the “dipole approximation” for subsurface scattering?

If you read papers about subsurface scattering, you'll frequently come across references to something called the "dipole approximation". This term seems to go back to the paper A Practical Model for Subsurface Light Transport by Henrik Wann Jensen et al, but this paper is pretty difficult to understand.

Can anyone explain in relatively simple terms what the dipole approximation is and how it's used in rendering subsurface scattering?

• It is indeed a complicated paper (still gives me nightmares from time to time!) but I tried to simplify the whole thing a bit. Let me know if you think the below answer should be adjusted in a way or another :) – cifz Aug 6 '15 at 8:29

The assumption underlying such model is the same as lots of other models for skin rendering; the subsurface scattering can be approximated as a diffusion phenomenon. This is good because in highly scattering media, the distribution of light loses dependency from the angle and tends to isotropy.

The dipole approximation is a formulation for the resolution of such diffusion problem in an analytical fashion.

Basically they start by approximating the BSSRDF as a multiple scattering and single scattering component. The multiple scattering is then defined as:

Where $F_t$ are Fresnel terms and $R$ is the diffusion profile expressed as function of the distance between the entry and exit point.

This $R$ is referred to as diffusion profile and they formulate this profile via a dipole approximation. The contribution of the incoming light ray is considered to be the one of two virtual sources: one negative beneath the surface and one positive above it (that's why dipole)

Here in the picture r is the $\|x_i - x_o\|$ above. The contribution of those light sources is dependent on various factors such as the distance of the light from the surface, scattering coefficient etc. (See below for a more detailed description of the formula itself).

This model only account for multiple scattering events, but that's fine enough for skin. It must be noticed though that for some translucent materials (e.g. smoke and marble) the single scattering is fundamental. That paper propose a single scattering formulation, but is expensive.

The diffusion profile is usually approximated for real-time application as a series of gaussian blurs (like in the seminal works of D'Eon et al. in GPU Gems 3 then used for the Jimenez's SSSSS) so to make it practical for real time scenarios. In this wonderful paper there are details on such approximation. A picture from that paper show actually how good is this formulation:

As side note the dipole approximation assumes that the material is semi-infinite, however this assumption doesn’t hold with thin slabs and multi-layered material such as the skin. Building on the dipole work, Donner and Jensen [2005] proposed the multi-pole approximation that accounts for the dipole problems. With this model instead of a single dipole, the authors use a set of them to describe the scattering phenomenon. In such formulation the reflectance and transmittance profiles can be obtained by summing up the contribution of the different dipoles involved

EDIT: I am putting here the answers to a couple of @NathanReed 's questions in the comment section:

Even with the diffusion profile approximation, the BSSRDF model still requires integrating over a radius of nearby points on the surface to gather incoming light, correct? How is that accomplished in, say, a path tracer? Do you have to build some data structure so you can sample points on the surface nearby a given point?

The BSSRDF approximation still need to be integrated over a certain area, yes.

In the paper linked they used a Montecarlo ray-tracer randomly sampling around a point with a density defined as:

$\sigma_{tr}e^{-\sigma_{tr}d}$

Where that sigma value is the effective extinction coefficient defined below ( it is dependent on scattering and absorbing coefficient, which are properties of the material) and d is the distance to the point we are evaluating. This density is this way defined because the diffusion term has exponential fall-off.

In [Jensen and Buhler 2002] they proposed an acceleration technique. One of the main concepts was to decouple the sampling from the evaluation of the diffusion term. This way they perform a hierarchical evaluation of the information computed during the sampling phase to cluster together distant samples when it comes to evaluating the diffusion. The implementation described in the paper uses an octree as structure. This technique, according to the paper, is order of magnitude faster than the full Monte Carlo integration.
Unfortunately I never got myself into an off-line implementation, so I can't help more than this.

In the real-time sum-of-Gaussians approximations the correct radius is implicitly set when defining the variance of the Gaussian blurs that need to be applied.

Why one positive and one negative light? Is the goal for them to cancel each other in some way?

Yes, the dipole source method (which dates way before Jensen's paper) is such defined to satisfy some boundary condition. Specifically the fluence must be zero at a certain extrapolated boundary that has a distance from the surface of $2AD$ where

Being $F_{dr}$ the fresnel reflectivity of the slab considered and that sigma value is the reduced extinction coefficient described below.

EDIT2: I have expanded (a tiny bit) some of the concepts in this answer in a blog post: http://bit.ly/1Q82rqT

For those who are not scared by lots of greek letters in a formula, here's an extract from my thesis where the reflectance profile is briefly described in each term:

• Thanks! I'm still wondering a couple of things. (1) Even with the diffusion profile approximation, the BSSRDF model still requires integrating over a radius of nearby points on the surface to gather incoming light, correct? How is that accomplished in, say, a path tracer? Do you have to build some data structure so you can sample points on the surface nearby a given point? And (2) Why one positive and one negative light? Is the goal for them to cancel each other in some way? – Nathan Reed Aug 6 '15 at 18:21
• 1) Indeed, what they propose in the paper with their monte carlo ray tracer is a stochastic sampling with a specific density based on distance and extinction coefficient. I guess you can dart-throw to find samples and use an appropriate acceptance probability based on the extinction coeff. and distance. (1/2) – cifz Aug 6 '15 at 19:37
• I know that Jensen published a hierarchical integration approach in 2002 that unfortunately I just read once and a while ago, so I just remember few bits. The core concept was to decouple the sampling from the diffusion approximation and cluster distant samples. IIRC they used an octree as hierarchical structure. I never got myself into an offline implementation so I am not of so much help on other details on this I am afraid. (2/2) – cifz Aug 6 '15 at 19:38
• 2) It is this way to satisfy some boundary conditions, you want to let the fluence become zero at a certain extrapolated boundary that has a specific distance from the medium. This distance is computed based on scattering coefficient and scattering albedo. – cifz Aug 6 '15 at 19:38
• @NathanReed Let me know if this has clarified something, otherwise I can try and expand the thoughts in this comment in the answer itself – cifz Aug 7 '15 at 5:38

To easily understand the 'dipole theory' we have first to understand from where it comes from the 'diffusion theory'.

And it comes from simulating light transport in participating media by solving the radiance transport equation (RTE).

The classical diffusion approximation solves the RTE by considering only a first-order spherical harmonic expansion of the radiance. If we further assume that the source function is a unit power isotropic point source in an infinite homogeneous medium we get to the classical diffusion Green's function (monopole).

When rendering translucent materials, it is convenient to re-formulate this problem in analogy to the local surface reflection integral. This results in an integral equation which computes the outgoing radiance, Lo, at position and direction as a convolution of the incident illumination, Li, and the BSSRDF, S, over all incident positions and directions.

Now with Green's diffusion function we don't account for boundary conditions imposed by the material surfaces. These conditions can be simply handled by placing a mirrored negative source outside the medium for every positive source inside the medium such that the fluence is zeroed at an extrapolated distance above the surface. This is the dipole approximation.

So without boundary conditions our fluence is expressed volumetrically within the medium. To render subsurface scattering we have to compute light leaving various points on the surface. For this to happen we have to compute diffusion profiles due to a dipole so to evaluate the directional derivative of the fluence in the direction of the surface normal.

Latest advancements on BSSRDF comes from shifting the initial volumetric and participating media approach with ad-hoc stuff that better takes into account surface boundary conditions.

As for the ..

BSSRDF model still requires integrating over a radius of nearby points on the surface to gather incoming light .. ?

Yes, we integrate incident light togheter with the BSSRDF over all incident positions and directions.

Now we can adopt here a brute force method or a dart-throwing approach with russian rouletting. But they are both naive approaches.

SSS came to production (Pixar Renderman) when an approximation method has been developed to use hierarchical light gathering using an octree data structure where to bake diffuse lighting to be available in clusters for SSS computation. With Renderman this was a natural approach because of REYES so that every micropoligon generated from REYES could be easily 'splatted' into a point and inserted into an octree.

A second wave of improvements relies on disk based importance sampling (Arnold) that is now the de-facto approach for many SSS implementations outthere. Generally we define a search volume (sphere), distribute samples on disk above the surface and probe along normal directions and orthogonal directions to find all hits inside the volume.