function shit = ellipseMC

f1 = @(x,y) x.^2/16 + y.^2 - 1;
fimplicit(f1)
hold on

Xrand = 8.*rand(10000,1) -4;
Yrand = 2.*rand(10000,1) -1;

Nin = 0;
Nout = 0;

plot(Xrand,Yrand,'.');

for ii = 1:10000
    if ((Xrand(ii)^2/16 + Yrand(ii)^2 <= 1))
        Nin = Nin + 1;
    else
        Nout = Nout + 1;
    end
end

disp(Nin/Nout);