% Compare convergence of random sampling with stratified sampling
% for Monte Carlo integration of the 3-dimensional unit hypersphere

% According to Eric Veach, "Robust Monte Carlo Methods for Light Transport
% Simulation", Ph.D. dissertation, Stanford University, December 1997,
% Section 2.6.1, stratified sampling is helpful only for small Ndim.
%     http://graphics.stanford.edu/papers/veach_thesis/
% Convergence rates for stratified sampling vs. Ndim were obtained by:
% Mitchell, D. P. [1996], "Consequences of stratiﬁed sampling in computer 
% graphics", SIG-GRAPH 96 Proceedings, Addison-Wesley, pp. 277–280.
%     O( N^(-(1+2/Ndim)/2) ) for integrands with bounded first derivative,
%     O( N^(-(1+1/Ndim)/2) ) for piecewise continuous integrands.
% The example here (Ndim=3 dimensions, integrand with infinite first
% derivative at one boundary) is consistent with the second result.


Ndim=3  % display number of dimensions for reference
correct = answer(Ndim); % exact volume, for comparison
p = correct/(2^Ndim) % expected fraction of points inside sphere
Ntrials = 100
for m = 1:5
   Npoints(m) = 8^(m+1);

   for trial = 1:Ntrials
      err(trial) = (hypersphere1(Ndim, Npoints(m)) - correct)/correct;
         % fractional error
   end
   err1(m) = sqrt(mean(err.^2));
   
   for trial = 1:Ntrials
      err(trial) = ( threesphere(Npoints(m)) - correct )/correct;
         % fractional error
   end
   err2(m) = sqrt(mean(err.^2));

end
figure(3)
hold off
loglog( Npoints, err1, 'ko', ...
        Npoints, err2, 'b^', ...
        Npoints, Npoints.^(-1/2), 'r-', ...
        Npoints, Npoints.^(-2/3), 'g-' );
legend( 'Random sampling', 'Stratified sampling', ...
   'N^{-1/2}', 'N^{-2/3}' )
title( 'Convergence: MC integration of the 3-sphere', 'fontsize', 20 )
xlabel( 'N, number of M-C samples', 'fontsize', 18 )
ylabel( 'RMS error', 'fontsize', 18 )