The current SampleHemisphere function does not generate points distributed uniformly on the hemisphere: 
Here's an efficient way to generate uniform samples on the unit sphere attributed to Marsaglia (1972).
You can flip the sample's position according to its dot product with a normal as before.
Here's an efficient way to generate uniform samples on the unit sphere attributed to Marsaglia (1972).
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 | internal v3
SampleSphere(random_series *Series)
{
v3 Result = {};
f32 X1, X2, SquaresSum, SquareRootTerm;
do
{
X1 = RandomBilateral(Series);
X2 = RandomBilateral(Series);
SquaresSum = Square(X1) + Square(X2);
} while (SquaresSum >= 1.0f);
SquareRootTerm = 2*SquareRoot(1 - SquaresSum);
Result.x = X1*SquareRootTerm;
Result.y = X2*SquareRootTerm;
Result.z = 1 - 2*SquaresSum;
return(Result);
}
|
You can flip the sample's position according to its dot product with a normal as before.
Here's a method which involves no loops, due to Archimedes (225 BCE). The idea is to randomly sample the enclosed cylinder then map that onto the sphere or hemisphere.
However, you probably don't actually want to sample the hemisphere uniformly; for light transport, you probably want to include the geometric cosine factor. The easiest way to do this is to randomly sample the unit disk, then map that onto a hemisphere sitting above it.
1 2 3 4 5 6 7 8 9 10 11 12 | internal v3
SampleHemisphere(random_series *Series)
{
v3 Result = {};
f32 z = RandomBilateral(Series);
f32 r = SquareRoot(1 - z*z);
f32 theta = RandomBilateral(Series) * 2 * PI;
Result.x = r * Cos(theta);
Result.y = r * Sin(theta);
Result.z = z;
return (Result);
}
|
However, you probably don't actually want to sample the hemisphere uniformly; for light transport, you probably want to include the geometric cosine factor. The easiest way to do this is to randomly sample the unit disk, then map that onto a hemisphere sitting above it.
1 2 3 4 5 6 7 8 9 10 11 12 13 | internal v3
CosineSampleHemisphere(random_series *Series)
{
v3 Result = {};
f32 r2 = RandomBilateral(Series);
f32 r = SquareRoot(r2);
f32 theta = RandomBilateral(Series) * 2 * PI;
f32 z = SquareRoot(1 - r2);
Result.x = r * Cos(theta);
Result.y = r * Sin(theta);
Result.z = z;
return (Result);
}
|