hist_samples = 1e6;
K = 12;
mu = K/2;
sigma2 = K/12;
% Should look like a standard normal for large K:
samples = (sum(rand(K, hist_samples), 1) - mu) / sqrt(sigma2);

clf; hold on;
[hn, hx] = hist(samples, 50);
bar(hx, hn, 1);
bar_width = hx(2) - hx(1);
xx = hx(1):bar_width/5:hx(end);
% Multiply probablity density by bar_width to get bin probability, and further
% multiply by hist_samples to get the expected number of histogram counts:
plot(xx, hist_samples*bar_width*exp(-0.5*xx.^2)/sqrt(2*pi), 'r');


% A load of plotting junk:
title(sprintf('approx randn with %d uniforms', K));
xlabel('x');
ylabel(sprintf('counts out of %g draws', hist_samples));
box on
set(gca,'XColor', 0.4*ones(1,3))
set(gca,'YColor', 0.4*ones(1,3))
set(gca,'TickDir', 'out')
set(gca,'XMinorTick', 'on')
figname = 'clt_demo';
set(gcf, 'PaperPosition', [0 0 6 3]); % Makes the plot "chunkier"
print([figname, '.eps'], '-depsc');
system(['epstopdf ' figname, '.eps']);