
%%
%%%%%%%%%%%% 1D random numbers
N = 500000; 
X = rand (1,N) ;
hist(X,5) ; % shows distribution of X, probability mass function

X = rand (1,N) ;
Z = (X + Y)/2; 
hist(Z,5) ; % pmf of Z is not flat!

%%
%%%%%%%%%%%% 2D random numbers
N = 5; 
X = rand (2,N) ;


%% 2D integration
N = 50; 
X = rand (2,N) ; % points
r = 0.1; % radius
C = [0.65 0.25] ; % centre of circle

% assume none of the points is inside
inside = zeros(1,N) ; % array of flags

for i = 1:N
    if dist(C,X(:,i)) < r
        inside(i) = 1; 
    end
end

numinside = sum (inside) ;
CircleAreaEstimate = numinside/N ;
piestimate = CircleAreaEstimate/(r*r) 

% plot samples and circle
clf;
axis([0,1,0,1]);
hold on;
viscircles(C, r, 'color', 'blue') ;
scatter(X(1,:), X(2,:), 'k', 'markerfacecolor', 'black');
for i=1:N
 if inside(i)
     scatter(X(1,i), X(2,i), 'k', 'markerfacecolor', 'red');
 end
end
hold off;

%% 
%%%%%%%%%%%%%%%%%% Repeat estimation of pi


%%
%%%%%%%%%%%%%%% arbitrary function
f = @(x) 5*x.^4 - 4*x.^2 + 2*x + 1.5 ;
xreg = linspace(0,1,10) ;
plot(xreg,f(xreg)) ;

N = 100 ;
X = rand(1,N);
I = mean(f(X)) 
