Test cosine weighted sampling and stratified sampling

testSampling.m

function testSampling(n)

samples = zeros(n, 3);
samples2 = zeros(n, 3);

for i=1:n
    xi1 = rand();
    xi2 = rand();
    samples(i, :) = hemi(xi1, xi2);
    samples2(i, :) = hemi_stratified(xi1, xi2);
end

figure; hold on;
plot3(samples(:,1), samples(:,3), samples(:,2), 'bx');
plot3(samples2(:,1), samples2(:,3), samples2(:,2), 'ro');
grid on;

end

function S = hemi(xi1, xi2)

SQRT_OF_ONE_THIRD = sqrt(1.0 / 3.0);
TWO_PI = 2.0 * pi;
up = sqrt(xi1);
over = sqrt(1 - up * up);
around = xi2 * TWO_PI;

normal = [0, 1, 0];

if abs(normal(1)) < SQRT_OF_ONE_THIRD
    directionNotNormal = [1, 0, 0];
elseif abs(normal(2)) < SQRT_OF_ONE_THIRD
    directionNotNormal = [0, 1, 0];
else
    directionNotNormal = [0, 0, 1];
end

perpendicular1 = cross(normal, directionNotNormal);
perpendicular1 = perpendicular1 / norm(perpendicular1);
perpendicular2 =            cross(normal, perpendicular1);

S = ( up * normal ) + ( cos(around) * over * perpendicular1 ) + ( sin(around) * over * perpendicular2 );

end

function S = hemi_stratified(xi1, xi2)

SQRT_OF_ONE_THIRD = sqrt(1.0 / 3.0);
TWO_PI = 2.0 * pi;
up = 1-xi1;
over = sqrt(xi1*(2.0-xi1));
around = xi2 * TWO_PI;

normal = [0, 1, 0];

if abs(normal(1)) < SQRT_OF_ONE_THIRD
    directionNotNormal = [1, 0, 0];
elseif abs(normal(2)) < SQRT_OF_ONE_THIRD
    directionNotNormal = [0, 1, 0];
else
    directionNotNormal = [0, 0, 1];
end

perpendicular1 = cross(normal, directionNotNormal);
perpendicular1 = perpendicular1 / norm(perpendicular1);
perpendicular2 =            cross(normal, perpendicular1);

S = ( up * normal ) + ( cos(around) * over * perpendicular1 ) + ( sin(around) * over * perpendicular2 );

end