function main()
%	montecarloarea(@myfunc1,-1,1,-1,1,10000)
	montecarloarea(@myfunc2,0,1,0,1,10000)
end

function f=myfunc1(x,y)

   f = x.^2 + y.^2 - 1;
end

function f=myfunc2(x,y)

   leftside = 4*(2*x-1).^4 + 8*(2*y-1).^8;
   riteside = 1 + 2*(2*y-1).^3 .* (3*x-2).^2;
   f = leftside - riteside;

end

function a=montecarloarea(func,xmin,xmax,ymin,ymax,n)
   x = xmin + (xmax-xmin).*rand(1,n);
   y = ymin + (ymax-ymin).*rand(1,n);
   f = func(x,y);
   negf = (f<0);
   count = norm(double(negf),1)  % Must be a better way of counting the 1s!
   a = (count/n)*(xmax-xmin)*(ymax-ymin);

   % make a picture
   x1 = x(f<0);
   y1 = y(f<0);
   hold off
   plot(x1,y1,'.b')
   axis equal
   hold on
   x1 = x(f>0);
   y1 = y(f>0);
   plot(x1,y1,'.r')
   axis equal

end