A Hexagonal Aperture Would Break Every Diffuse Surface in My Renderer

I have experimented with a lot of features on my raytracing journey. I played with reflections, implemented shadows, and lost braincells getting scenes to render. Finally, I implemented defocus blur.

The image reminded me of the bokeh you see in photos. Sometimes they aren't always round.

So I tried it with what I already had: two passes of the same scene, one circular and one hexagonal.

At that point I felt like I could turn anything into bokeh. I could have hearts, stars, concentric circles. Why not push it further?

Then I noticed something. Every diffuse surface in my renderer would suffer.

Here is the primitive that both my lens sampler and diffuse surface lean on. It rejection-samples a square until a point lands inside the unit disk, and that test is the only reason the shape comes out round.

inline glm::vec2 SampleUnitDisk(Rng& rng) {
  while(true){
    float x = rng.Range(-1.0f, 1.0f);
    float y = rng.Range(-1.0f, 1.0f);
    if (x * x + y * y <= 1.0f) return glm::vec2(x, y);
  }
}

And here is what the hexagonal change looked like:

// Regular hexagon inscribed in the unit circle: corners at radius 1, flat top
// and bottom, flat edges at the apothem sqrt(3)/2 ~= 0.866.
constexpr float kApothem = 0.8660254f;

inline glm::vec2 SampleUnitDisk(Rng& rng) {
  while (true) {
    float x = rng.Range(-1.0f, 1.0f);
    float y = rng.Range(-1.0f, 1.0f);
    // WAS: if (x * x + y * y <= 1.0f) return glm::vec2(x, y);
    if (std::fabs(y) <= kApothem &&
        std::fabs(kApothem * x + 0.5f * y) <= kApothem &&
        std::fabs(kApothem * x - 0.5f * y) <= kApothem) {
      return glm::vec2(x, y);
    }
  }
}

The two separate call sites:

// The camera: a physical aperture.
const glm::vec2 lens_uv = SampleUnitDisk(rng) * cfg.aperture_radius;

// The diffuse BSDF wrapper
inline glm::vec3 CosineWeightedHemiSphereSurface(Rng& rng,
                                                 const glm::vec3& normal,
                                                 const glm::vec3& tangent,
                                                 const glm::vec3& binormal) {
  const glm::vec2 d = SampleUnitDisk(rng);
  const float z = std::sqrt(std::max(0.0f, 1.0f - d.x * d.x - d.y * d.y));
  return d.x * tangent + d.y * binormal + z * normal;
}

Why do my BSDF and my lens share a disk sampler?

Fair question. If testing, or more fittingly, playing with my lens can break my diffuse surfaces, why the coupling?

The BSDF side is Malley's method. Sample a point inside a unit disk, then project it along the normal onto a hemisphere originating from the ray's intersection with the surface. That is how you implement a diffuse surface.

The lens side never heard of Malley. Randomizing around the origin allows the ray to focus on a particular plane. Near and far objects blur out, while whatever sits on that plane stays in focus. Guess how we randomize around the origin. Exactly!! We sample within a unit disk.

So both sides have to sample a unit disk, and voila; one function, two callers. Apparently PBRT does it the same way, so I will take that as great minds thinking alike.

I only caught it by reading the code

So how did I catch it? I read the code, because I wrote it. Nonetheless, I cannot rely on that catch. I could easily have missed it, someone opening the file a year from now could miss it too, and neither of us would know the diffuse was off until much later.

Looking at the render does not help either. I can make the hexagon, sure, but the diffuse barely moves.

Do you notice the difference? Do you see the shadows shift? Yeah, neither do I. Without the two frames side by side I would have walked straight past this.

So I said to myself, "What can I actually check for?"

Check that it is still a circle? I can only answer that by looking at what the sampler draws rather than the lines that produce it.

Any sampler is correct as long as the Pdf() tells the truth

So start with the shape mine draws, and why it works at all. It works because solid angle projects onto the disk with a factor of cosine, $dA = \cos\theta \, d\omega$ , so a uniform disk sample lifted to the hemisphere has density

For a Lambertian surface, whose BRDF is $\rho/\pi$ , the estimator then cancels completely:

Each bounce multiplies output by exactly the albedo, ρ (the fraction of light a surface sends back, [0,1]).

Read that derivation again, though, and notice what it never asks for. It never asks for a disk.

Which means I can build the same lobe without one. Here is the road Ray Tracing in One Weekend takes: stand at the tip of the normal, push out by a random unit vector, and normalise where you land.

const glm::vec3 dir = glm::normalize(normal + RandomUnitVector(rng));

That one draws from a ball. It never touches a disk, and it arrives at the same cosine lobe anyway. Two roads, one distribution and I can compare without knowing in advance where the distribution is supposed to be. Get it memorized.

The lift, and what the hexagon does to it

Malley's method is one picture. Put the shape underneath the hemisphere like a floor plan, and let every point on it rise straight up to the dome. The height it lands at is cos θ. The cloud of directions you end up with is the lobe: the fat dome.

Look at the hexagon in that diagram and notice how little moved. Every point still lands above the horizon, the centre still maps to the normal, and it is still a fat lobe crowding around it.

What moves is the rim, and only the rim. A hexagon's corners reach r = 1, but its edges cut in to the apothem (as close as an edge ever gets to the centre) = √3/2 ≈ 0.866, and √3/2 lifts to cos θ = 0.5. It's SOH CAH TOA again. So the grazing directions survive in six chunks and are gone the rest of the way round.

Which is enough to price the damage before anything runs, because the lift hands every question about the lobe back to the floor plan. Malley over a region $R$ of area $A$ gives density $\cos\theta / A$ , so every moment is an integral over the shape:

Read $\mathbb{E}[\,\cdot\,]$ as "the average of". At $k = 2$ we lose the square root entirely:

so the second moment is the floor plan's mean squared radius. At $k = 1$ , one polar integral per edge, and for a hexagon it closes:

Nothing is specific to hexagons. Hand it a square, a heart, a five-pointed star, and the same two integrals price that sampler too. And the other road has no floor plan at all, which is the point of it: built out of SampleUnitBall, it sits on the disk's row, 2/3 and 1/2, exactly, whatever shape the lens becomes.

So we currently have a few predictions in before a single test run: 0.744 against 0.667 with a gap of 0.077.

Get that memorized, also.

So I built the bug on purpose

I reverted to a hexagon in SampleUnitDisk and ran my tests.

The furnace passed, as expected. All 10 render tests stayed green, the 6 furnace tests among them. The furnace check essentially checks the radiance, and that was virtually no different.

The two below it check shape (for a Lambertian the second is π times the first on average; the two tests draw separate samples) and both failed. In line with our prior prediction. Remember the 0.744 number. Yeah, we are in business now.

The test that caught it

Let's dive deeper into that test. Before making the change I had already written a test that builds the cosine lobe that way, and it draws from SampleUnitBall, not SampleUnitDisk, so it cannot move when the aperture does. Two constructions of one distribution have to agree on every moment, which makes them a differential test: no reference renderer, no known answer needed. Over 2 million samples each:

We can see from the table above that we were right to 3 significant figures, and the gap that is left is within sampling noise. In other words, we were still right. The earlier predictions for both E[cos] and E[cos²] match up.

What I do differently now

• Share on mechanism, never on coincidence. I will still share primitives. The ideal move is creating new primitives for specific cases. I already knew that. Utilizing a proper separation of concerns.
• Test the shape, not just the total. The furnace passed, and it measures totals. I still needed to test the shape.
• Build the bug on purpose. You just might learn something.

The shapes I wanted in the first place

The lesson was never "don't shape the aperture." It was "give the lens its own sampler first." So that is what I did, and then I spent an evening on the part that has no excuse.

Look at all those images. "chef kiss." This was really all I wanted and it could have broken such a key section of my ray tracer. Sometimes that is the price of having one.

From a path tracer I'm building in C++. Next: getting total internal reflection wrong in four different ways before getting it right.