0
$\begingroup$

I've been looking at Monte Carlo with PBRT. Currently I'm aiming to finish up the basic MC so that I can move on to multiple importance sampling. For diffuse surfaces I'm calculating the random sample direction with cosine weighting i.e. normalize(normal + uniformSphere()). PBRT states that I should divide each sample by the PDF, and as I'm using cosine weighting, the PDF should be cos(theta) or simply dot(normal, direction). The major issue I'm observing is that this value is always less than 1, so every sample ends up brighter than it should be. The test case I'm using has no actual light sources, only ambient light, so every sample produces the same constant result, however the image now has a large amount of noise where none existed before, because the ambient light can be sampled from many different directions which the PDF weights differently.

How do I make sure that weighting the samples doesn't impact their final brightness? I've seen various sources on this but none of them seem to have any kind of renormalisation which you would normally see for a weighted average.

$\endgroup$

2 Answers 2

3
$\begingroup$

There are several misunderstandings:

  • normalize(normal + uniformSphere()): I don't think this is the correct way to implement cosine-weighted hemisphere sampling. Check this out: Sampling the hemisphere and 13.6.3 Cosine-Weighted Hemisphere Sampling (this is more involved). The method you implemented by normalizing the sum of normal and uniformly sampled vector is neither uniform hemisphere sampling nor cosine-weighted sampling.
  • PDF is not cos(theta). PDF is cos(theta) / pi. If you integrate the PDF over the hemisphere then you will find that the latter one results in 1 (valid). Don't forget to add the term $\sin \theta$ to account for the measure conversion between solid angle and polar coords.

Getting brighter results is caused by not using unbiased sampling method and PDF pairs. Let's break it down (a little bit of math, but I will try to make it easy to understand):

MC integration is aimed at approximating the following integral in an unbiased way:

$$ \int_{\Omega} f(x)d\mu(x)\tag{1} $$

That is, for the following estimator (let's consider a one-sample case, N = 1): $$ \hat{I} = \underbrace{f(x)}_{\text{evaluation}} / \underbrace{p(x)}_{\text{some PDF you choose to divide}}\tag{2} $$ If we take its expectation: $$ \mathbb{E}(\hat{I}) = \int_{S}\frac{f(x)}{p(x)}p_{act}(x)d\mu(x)\tag{3} $$

We want $p(x)$ (the PDF you choose as the denominator) and $p_{act}(x)$ (the actual PDF of the samples, determined by the sampling method you use) can cancel each other out, so that the result will be exactly $(1)$. Now that we know this, for the specific case of yours:

If you choose $p(x)$ as $\cos$ or even $\cos / \pi$: since the reflected direction given by your code is not drawn from the cosine-weighted PDF, $p_{act}$ won't cancel itself out with $p$ in the denominator, which results in bias (hence, brighter). Another extreme example is that if you choose $\cos /\pi$ as $p$ and use uniform hemisphere sampling, $p$ and $p_{act}$ won't cancel each other out, either. So, to get the desired output, always choose sampling method ($p_{act}$) and the denominator $p(x)$ wisely.

What you should do here, is to:

  1. Implement the local reflection direction sampling, according to Cosine-weighted hemisphere.
  2. Rotate the local vector to the global frame: since you know the normal, and your cosine-weighted sample is drawn w.r.t to vector (0, 0, 1): calculate the rotation from (0, 0, 1) to your normal and apply the rotation to your local vector.
  3. Devide $f(x)$ by $\cos\theta / \pi$, the $\cos\theta$ term should cancel the foreshortening term in $f(x)$, so you will end up with a simple path throughput: $k_d$.
$\endgroup$
7
  • $\begingroup$ The thing that confuses me is that I don't see where the P act term comes from. In this case, f(x) just shoots a ray in that direction, so for my test scene f(x) = 0.1. So regardless of the probability of the ray in any given direction, the resulting samples aren't actually numerically changed, whereas dividing them by the PDF does change the result. $\endgroup$
    – Puppy
    Commented Jul 22 at 18:41
  • $\begingroup$ $p_{act}$ is given by the actual distribution of the sample. The example is given in the post: if you sample the hemisphere by uniform sampling, then $p_{act}$ will be the PDF of uniform hemisphere distribution, regardless what you choose as the denominator. $\endgroup$
    – Enigmatisms
    Commented Jul 23 at 1:19
  • $\begingroup$ So, if you choose uniform hemisphere sampling to generate samples, while divide $f(x)$ by its PDF, you will have unbiased results. The same holds for cosine-weighted sampling. To understand why this is the case, always take the expectation of the stochastic estimator, as in Equation (3). $\endgroup$
    – Enigmatisms
    Commented Jul 23 at 1:21
  • $\begingroup$ You will get $\int f(x) p_{act}(x) d\mu(x)$ as expectation, if you don't account for the actual distribution of samples in your estimator. Yet what we desire, is $\int f(x) d\mu(x)$. $\endgroup$
    – Enigmatisms
    Commented Jul 23 at 1:23
  • $\begingroup$ I don't understand how that works. If I choose $f(x) = 0.1$ then all distributions produce the same result, so there's no guarantee that the actual distribution impacts the final result. Even in a more realistic case where I might get different results, it won't change e.g. the range. If $f(x)$ has the range 0-1 and $cos(\theta)$ has a maximum of 1, then $f(x)p_{act}(x)$ would have the range 0-$\frac{1}{\pi}$. But no matter how I distribute my samples, they still have the range 0-1. So changing the distribution doesn't seem to have the same properties as multiplying the sample value. $\endgroup$
    – Puppy
    Commented Jul 23 at 18:38
0
$\begingroup$

It seems like the confusion comes from some assumptions that go unstated. The division by the PDF is given as a general rule for all Monte Carlo, but this doesn't actually seem to be true. In practice, PDFs which cancel out with part of f are used, not any arbitrary PDF, even for functions where this should work like f(x) = 0.1. PBRT is a little unhelpful on this front because first the division by PDF is given, then later on, the parts of f that it should cancel with. So it seems like the correct way of describing it should be that given that f(x) = g(x)p(x), then sampling g(x) by p(x) and cancelling out p(x) is valid. Otherwise, it is not valid.

$\endgroup$
1
  • $\begingroup$ It is true, $\int_{\Omega}f(x) \,dx \approx \frac{1}{N}\sum_{k=1}^N \frac{f(x_k)}{p(x_k)}$ is valid as long as $x_1, \ldots, x_N$ are independently and identically distributed with density $p$. It doesn't matter whether $f$ was formed as $f(x) = p(x) g(x)$ originally. You can rewrite any function like this anyways. $\endgroup$
    – lightxbulb
    Commented Aug 1 at 17:18

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.