function ezisosurf(f,ranges,fval,m)
% ezisosurf(f,[x1 x2 y1 y2 z1 z2],[f1 ... fn])
% Evaluate the function f on grid for x in [x1,x2], y in [y1,y2], z in [z1,z2].
% Then draw contour surfaces for f=f1,...,f=fn

if nargin==3
  m = 21;
end
ih = ishold;
xval = linspace(ranges(1),ranges(2),m);
yval = linspace(ranges(3),ranges(4),m);
zval = linspace(ranges(5),ranges(6),m);
[X,Y,Z] = meshgrid(xval,yval,zval); % values of x,y,z on grid
ffun = inline(vectorize(f),'x','y','z');     % turn symbolic expression into function
fg = ffun(X,Y,Z);                            % values of f on the grid
for val = fval                               % plot level surfaces for these values of f
  isosurf(X,Y,Z,fg,val); hold on;
end
if ~ih
  hold off
end

