implementation:

1. Microfacet material
2. Multiple importance sampling
3. Malley's Method
4. reflection and refraction

Taichi Course Season 1 final project (homework).

Environment:

[Taichi] version 0.8.4, llvm 10.0.0, commit 895881b5, win, python 3.7.10
[Taichi] Starting on arch=cuda

Reference: 《Physically based rendering:from theory to implementation》

path tracer based on 《PBRT》

一.introduction to sampling theory

1. what is sampling?

impulse train：

sampling process corresponds to multiplying the function by a “impulse train” function, an infinite sum of equally spaced delta functions.

《PBRT》A digital image is represented as a set of pixel values, typically aligned on a rectangular grid. When a digital image is displayed on a physical device, these values are used to determine the spectral power emitted by pixels on the display.

《PBRT》the pixels that constitute an image are point samples of the image function at discrete points on the image plane.

there is no “area” associated with a pixel.

when sampling the film signal

pos = camera_pos
ray_dir = ti.Vector([
            (2 * fov * (u) / resolution[1] - fov * resolution[0] / resolution[1] - 1e-5),
            2 * fov * (v) / resolution[1] - fov - 1e-5, -1.0
        ]).normalized()

then we need anti-aliazing

pos = camera_pos
ray_dir = ti.Vector([
            (2 * fov * (u + ti.random()) / resolution[1] - fov * resolution[0] / resolution[1] - 1e-5),
            2 * fov * (v + ti.random()) / resolution[1] - fov - 1e-5, -1.0
        ]).normalized()

二.sampling

Preview (CDF sampling technique)

There are many techniques for generating random variates from a specified probability distribution such as the normal, exponential, or gamma distribution. However, one technique stands out because of its generality and simplicity: the inverse CDF sampling technique.

1. Uniformly Sampling a Hemisphere (multidimensional sampling technique)

a uniform distribution means that the density function is a constant, so we know that p(x) = c

so p(ω) = 1/2*pi

then p(θ, φ) = sinθ/2*pi

Notice that the density function for φ itself is uniform

then use the 1D inversion technique to sample each of these PDFs in turn

2. sample area light

def sample_area_light(hit_pos, pos_normal):
    # sampling inside the light area
    x = ti.random() * light_x_range + light_x_min_pos
    z = ti.random() * light_z_range + light_z_min_pos
    on_light_pos = ti.Vector([x, light_y_pos, z])
    return (on_light_pos - hit_pos).normalized()

3. introduction to importance sampling

why we need importance sampling?

the Monte Carlo estimator converges more quickly if the samples are taken from a distribution p(x) that is similar to the function f(x) in the integrand.

《PBRT》：We will not provide a rigorous proof of this fact but will instead present an informal and intuitive argument.

then we try to analyze the importance sampling method

we have three terms

• BRDF
• incident radiance ( infeasible )
• cosine term

4. cosine-weighted sampling

Malley's method

So, We could compute the marginal and conditional densities as before, but instead we can use a technique known as Malley’s method to generate these cosine-weighted points.

1. cosine term

2. 2D Sampling with Multidimensional Transformations

   • (1) sampling a unit disk (Concentric Mapping)

   • (2) project up to the unit hemisphere (cosine-weighted hemisphere sampling)

(1) sampling a unit disk

(2) projection

To complete the (r,φ)=(sinθ,φ)⇒(θ,φ) transformation, we need the determinant of the Jacobian

Why

5. multiple importance sampling

BDPT only：

BDPT + MIS：

Why we need MIS?

• balance heuristic

• power heuristic (Veach determined empirically that β=2 is a good value.)

//Compute heuristic
def mis_power_heuristic(pf, pg):
    # Assume 1 sample for each distribution
    f = pf ** 2
    g = pg ** 2
    return f / (f + g)
    # return 1

//combine
@ti.func
def sample_light_and_cosineWeighted(hit_pos, hit_normal):
    cosine_by_pdf = ti.Vector([0.0, 0.0, 0.0])

    light_pdf, cosineWeighted_pdf = 0.0, 0.0

    # sample area light => dir, light_pdf; then dir => lambertian_pdf; then mis
    light_dir = sample_area_light(hit_pos, hit_normal)
    if light_dir.dot(hit_normal) > 0:
        light_pdf = compute_area_light_pdf(hit_pos, light_dir)
        cosineWeighted_pdf = compute_cosineWeighted_pdf(hit_normal, light_dir)
        if light_pdf > 0 and cosineWeighted_pdf > 0:
            l_visible = visible_to_light(hit_pos, light_dir)
            if l_visible:
                heuristic = compute_heuristic(light_pdf, cosineWeighted_pdf)
                DoN = dot_or_zero(light_dir, hit_normal)
                cosine_by_pdf += heuristic * DoN / light_pdf

    # sample cosine weighted => dir, lambertian_pdf; then dir => light_pdf; then mis
    cosineWeighted_dir = cosine_weighted_sampling(hit_normal)
    cosineWeighted_pdf = compute_cosineWeighted_pdf(hit_normal, cosineWeighted_dir)
    light_pdf = compute_area_light_pdf(hit_pos, cosineWeighted_dir)
    if visible_to_light(hit_pos, cosineWeighted_dir):
        heuristic = compute_heuristic(cosineWeighted_pdf, light_pdf)
        DoN = dot_or_zero(cosineWeighted_dir, hit_normal)
        cosine_by_pdf += heuristic * DoN / cosineWeighted_pdf

    # direct_li = mis_weight * cosine / pdf
    return cosine_by_pdf