Project 3-2: Pathtracer 2

CS 184: Computer Graphics and Imaging, Spring 2019

Andrew Campbell, cs184-adu

Overview

In this project, I extended my pathtracer to support more complicated materials, including glass and mirrors, and added environment lighting and depth of field effects. Whereas the first project was about setting up the core pathtracing infrastructure, this project was about modeling complex surfaces to produce truly striking results.

The main ideas in this project are conceptually simple; the main implementation challenge was in efficient sampling techniques. This project was at times frustrating because the calculations are complex, making debugging difficult. The results, however, are stunning and highly rewarding. In the process, I gained insight into the fascinating field of material modeling and gained a greater appreciation for abilities of pathtracing.

Part 1: Mirror and Glass

The main idea here is to implement the proper BSDF (bidirectional scattering distribution function). The BSDF is a generalization of the BRDF, as it additionally includes a BTDF, or bidirectional transmittance distribution function. We used a diffuse (constant) BSDF in the last project; now we support more complex ones.

Mirror material

We support mirror surfaces by implementing MirrorBSDF::sample_f(), which returns a sample of the surface BSDF. Given outgoing direction wo, we get wi by simply reflecting wo about the surface normal. This is trivial in the object coordinate space that BSDF calculations occur in; since the $z$ axis lies along the normal, we need only negate the $x$ and $y$ coordinates.

The pdf is just $1.0$ since it is a delta BSDF; there is no randomness involved. We return $\rho$, the reflectance, as the sampled value. Note that we scale the reflectance by the reciprocal of $\cos{\omega_i}$ so as to cancel out the cosine term used by the parent function at_least_one_bounce_radiance(); a perfect mirror does not produce Lambertian falloff.

Refraction

To implement materials like glass, we first need to support refraction. We implement this in BSDF::refract(), which calculates $\omega_i$ given $\omega_o$ in the object reference frame. Using Snell’s law, it can be shown that the correct equations are

$$\begin{align*} \omega_i.x &= -\eta \omega_o.x \\ \omega_i.y &= -\eta \omega_o.y \\ \omega_i.z &= \mp \sqrt{1 - \eta^2 (1-\omega_o.z^2)} \end{align*}$$

where $\eta$ is the ratio of the old index of refraction to new index of refraction; we compute $\eta$ in the code according to whether the ray is entering or exiting the non-air material. The $\mp$ indicates that we take the opposite sign of $\omega_o.z$. Note that if the radicand is negative, this is called total internal reflection.

Glass material

Using refraction, we can now implement GlassBSDF::sample_f(). Both reflection and refraction occur at an intersection with a glass material - we use a simple trick called Schlick’s approximation to give a coin-flip probability of either reflecting or refracting. By Schlick’s approximation, the ratio of reflection energy to refraction energy is

$$\begin{align} R(\theta) &= R_0 + (1 - R_0)(1 - \cos \theta)^5 \\ R_0 &= \left(\frac{n_1-n_2}{n_1+n_2}\right)^2 \end{align}$$

where $\theta$ is the angle between the incident light and the normal, and $n_1,\,n_2$ are the indices of refraction of the two media. We assume $n_1=1$ since one of the media is always air.

We do the same as a mirror BSDF if total internal reflection occurs or coinflip(R); otherwise we refract wo to *wi. For both cases, we weight the output Spectrum and pdf by the appropriate Schlick coefficient $R$ or $1-R$.

This part of the project was fairly straightforward and I didn’t run into any trouble.

Below we show a sequence of renderings of the scene CBspheres.dae with different max_ray_depth values. The left sphere has a mirror surface and the right sphere has a glass surface.

Beginning with max_ray_depth=0, only the light source is visible since all other light involves at least one bounce.

With max_ray_depth=1 (left), we see the effect of direct lighting. The spheres are dark because both reflection and refraction involve more than once bounce. With max_ray_depth=2 (right), the mirror ball appears because two bounces are needed for light to strike the ball and then the box.

With max_ray_depth=3 (left), the effect of refraction in the glass ball appears, since rays can enter the sphere, exit the sphere, and strike the box. The mirror reflection is 1 bounce behind, so to speak, so the glass ball reflection is dark. With max_ray_depth=4 (right), we see the caustic from light passing through the glass sphere since 4 bounces are needed (floor, into sphere, out of sphere, and box).

With max_ray_depth=5 (left), light can reflect off of the mirror ball and enter the glass sphere from the side, producing a small light patch on the blue wall. With max_ray_depth=100 (right), the multi-bounce effects have converged and the only noticeable difference is a slightly brighter scene due to the accumulation of indirect light.

Here is a high-res render of CBlucy.png that took about 3 hours.

---

Part 2: Microfacet Material

We now add support for certain microfacet models. A microfacet material is one whose surface is composed of many very tiny facets, each of which is a perfect specular reflector. In this project, we support isotropic rough conductors that reflect (not refract).

Microfacet BRDF

We implement the microfacet BRDF in MicrofacetBSDF::f(). This function is given by

$$f=\frac{F(\omega_i) * G(\omega_o, \omega_i) * D(h)}{4 * (n\cdot\omega_o) * (n\cdot\omega_i)}$$

$n$ is the surface normal, which is $(0,0,1)$ in object coordinates. $h = \frac{\omega_o + \omega_i}{\| \omega_o + \omega_i \|}$ is the half-vector. $D$ is the normal distribution function, which we implement in MicrofacetBSDF::D(); we use the Beckmann distribution given by

$$D(h)=\frac{e^{-\frac{\tan^2\theta_h}{\alpha^2}}}{\pi\alpha^2\cos^4\theta_h}$$

where $\alpha$ is the roughness of the macro surface (the smaller, the smoother) and $\theta_h$ is the angle between $n$ and $h$.

Lastly, $F$ is the Fresnel term implemented in MicrofacetBSDF::F(). It returns a Spectrum according to the approximation that each channel has a fixed wavelength:

$$F = \frac{R_s + R_p}{2} \\ R_s=\frac{(\eta^2+k^2)-2\eta\cos\theta_i+\cos^2\theta_i}{(\eta^2+k^2)+2\eta\cos\theta_i+\cos^2\theta_i} \\ R_p=\frac{(\eta^2+k^2)\cos^2\theta_i-2\eta\cos\theta_i+1}{(\eta^2+k^2)\cos^2\theta_i+2\eta\cos\theta_i+1}$$

where $\eta$ and $k$ represent indices of refraction for conducting materials; both are Spectrums encoding values at 614nm, 549nm, and 466nm. This website has a collection of $\eta$ and $k$ values for various materials at a given wavelength.

Importance Sampling

Cosine hemisphere sampling works well for a diffuse BRDF, but is highly inefficient for the Beckmann NDF due to its strongly preferential directions. Accordingly, we implement a form of importance sampling in MicrofacetBSDF::sample_f().

The implementation idea is to sample a $\theta_h$ and $\phi_h$ from some pdfs $p_\theta$ and $p_\phi$, respectively; combine them to form the sampled microfacet normal $h$; and reflect $\omega_o$ about $h$ to get the sampled $\omega_i$.

The pdfs we use resemble the Beckmann NDF:

$$p_\theta(\theta_h)=\frac{2\sin\theta_h}{\alpha^2\cos^3\theta_h}e^{-\tan^2\theta_h/\alpha^2} \\ p_\phi(\phi_h)=\frac{1}{2\pi}$$

We use the inversion method to sample them. We then recover $h$ by converting the $\theta$ and $\phi$ into $x, y, z$ coordinates. $\omega_i$ is the reflection of $\omega_o$ about $h$:

$$\omega_i = -\omega_o + 2 (\omega_o \cdot h) h$$

The pdf value is also calculated appropriately (derivation is tricky and omitted here).

We then return the BRDF evaluated at the sampled values.

I had some numerical issues in this part of the project, which I resolved by ensuring that all terms in the denominator of the BSDF are positive.

To illustrate the improvement gained from importance sampling, we compare two images of the scene CBbunny_microfacet_cu.dae rendered using cosine hemisphere sampling (left) and importance sampling (right). Both images were created with 64 samples per pixel, 1 sample per light, and 7 bounces.

Cosine hemisphere sampling gives a noisier result because it samples uniformly while importance sampling distributes samples where light is expected to have greater influence, thus producing better results for a given number of samples.

To illustrate the effect of $\alpha$ on the microfacet surface, we compare four renderings of CBdragon_microfacet_au.dae with different $\alpha$ values. All images were rendered with 1024 samples per pixel, 4 samples per light, and 7 bounces.

From left to right, top to bottom, we have $\alpha=0.005, 0.05, 0.25, 0.5$. The $\alpha$ value controls the roughness of the material, as given by the spread of the distribution of normals. The lower the value, the smoother the surface - resulting in a glossier appearance.

We can modify the $\eta$ and $k$ values to create any desired conductor material. Below, we attempt to model a bunny made out of mercury. The corresponding values (as given by the website above) are

Property	$614 \text{ nm}$	$549\text{ nm}$	$466\text{ nm}$
$\eta$	2.0733	1.1680	1.4612
$k$	5.3383	4.0572	4.5190

which results in

---

Part 3: Environment Light

In this part, we introduce environment light, which is an infinitely far away light source supplying radiance from all directions. The incoming light from each direction is defined in a texture map .exr file, which we parametrize by $\phi$ and $\theta$.

We do this in order to simulate more realistic lighting conditions, like outside settings on Earth where the entire sky is lit up from an extremely distant light source.

Environment Map

We first implement EnvironmentLight::sample_dir() to sample a Spectrum from an environment map at a particular direction. We simply use the helper functions to convert the direction vector into $\theta$ and $\phi$ coordinates, then $x$ and $y$ coordinates that can be used to index into the envMap. We then perform bilinear interpolation to sample a Spectrum at that point.

We first implement EnvironmentLight::sample_L() with uniform sampling, by getting a random direction on the sphere (i.e. pdf of $\frac{1}{4\pi}$) and converting it into environment map coordinates to be sampled with bilinear interpolation.

Importance sampling

We can improve the efficiency of light sampling given that most of energy is concentrated in the directions toward bright light sources.

The implementation idea is to convert the environment map into a 2D probability density function $p$ over $[0,w] \times [0,h]$. We represent this function as a 2D array that piecewise constant over the rectangle $[i,i+1] \times [j,j+1]$.

Without derivation, the pdf we use is

$$p(x,y) = \frac{E[\lfloor x \rfloor,\lfloor y \rfloor] \sin (\pi \lfloor y \rfloor/h)}{\sum_{i',j'} E[i',j'] \sin (\pi j'/h)}$$

where $E[i,j]$ represents an index into the environment map. Our first implementation step is to enforce

$$\int_0^w \int_0^h p(x,y) dx dy = \sum_{i=0}^{w-1} \sum_{j=0}^{h-1}p(i,j) \int_{j}^{j+1} \int_{i }^{i+1} dx dy = \sum_{i=0}^{w-1} \sum_{j=0}^{h-1} p(i,j) = 1$$

by properly normalizing by the sum in EnvironmentLight::init().

The sampling strategy is as follows: sample a row of the environment map using the marginal distribution $p(y)$; sample a pixel within that row using the conditional distribution $p(x \mid y)$; then convert that $(x,y)$ to a direction vector and return the appropriate radiance and pdf values.

We compute the cumulative marginal distribution in a double for loop as

$$F(j) = \sum_{j'=0}^j \sum_{i=0}^{w-1} p(i,j')$$

and store it in marginal_y. Then we compute the cumulative conditional distribution function as a double for loop as

$$F(i \mid j) = P(x < i | y = j) = \int_0^i p(x|y=j) dx = \int_0^i \frac{p(x,j)}{p(j)} dx = \sum_{i'=0}^{i-1} \frac{p(i',j)}{p(j)}$$

and store it in row-major form in conds_y.

We now update EnvironmentLight::sample_L() to perform importance sampling as follows: first, generate a uniform 2D sample on $[0,1]^2$. We get a row by finding the first j for which marginal_y[j] is greater than the first sample and we get a column by finding the first i for which conds_y[envMap->w * row + i] is greater than the second sample.

We convert the row and column into a direction ($\omega_i$) and then $x$ and $y$ indices with the same helper functions as earlier; and return the environment map value at $(x,y)$ using bilinear interpolation. We also make sure to compute the pdf (derivation tricky, omitted here).

For results, we will work with the following .exr file (shown here as converted to jpg):

The probability debug (visualization of the computed marginal and conditional distributions) is shown below.

I had a lot of trouble with debugging this part of the project, mostly due to an incorrect probability calculation. I had to take great care to convert the cumulative marginal distribution into the form expected by the cumulative conditional distribution equation.

The improvement importance sampling offers is illustrated below; both images are of bunny_unlit.dae with the environment map shown above with 4 samples per pixel and 64 samples per light.

The noise levels are noticeably lower in the importance-sampled render on the right, most notably in the bunny's shadow.

Doing the same comparison for bunny_microfacet_cu_unlit.dae, we have

Again, we see fewer speckles in the image produced by importance sampling.

Here is an additional rendering of bunny_microfacet_cu_unlit.dae in another environment map just for fun.

---

Part 4: Depth of Field

Up until now, our pathtracer used a pin-hole camera. All light entered through a virtual point, so everything was in focus. Real cameras, however, have finite apertures so that objects are in focus only if they’re in a plane focalDistance from the lens. We can use a thin lens approximation, i.e. ignore any effects due to the thickness of the lens, which is assumed to be negligible.

We use the following model to guide the implementation:

In this task, we generate the blue ray in the figure above. The implementation overview is to generate the red ray using the same technique as in the previous project; uniformly sample the disk representing the thin lens; calculate the focus point by using the focalDistance value and the red ray direction; generate a ray from the sampled disk point to the focus point; then properly perform camera-to-world conversion this ray.

The following renders show CBdragon.dae at various depths for a fixed lensRadius=0.01. From right to left, top to bottom, we show focalDistance=0.7, 0.9, 1.2, 2.5.

As we increase the focalDistance, the focus plane shifts from the tip of the dragon’s nose at 0.7 all the way to back corner at 2.5.

The following renders show CBlucy.dae at various apertures for a fixed focalDistance=1.0. From right to left, top to bottom, we show lensRadius=0.0, 0.2, 0.4, 0.6.

The entire scene is in focus at lensRadius=0.0, as this is just the pinhole camera case. As we increase the lensRadius, areas not at the focalDistance of 1.0 become increasingly out-of-focus.

Some additional renders of CBspheres_microfacet_al_ag.dae at different focalDistance values are shown for fun.

---