points=[-1 -1; 3 0; 0 4; 1 1; 2 5; 4 2]';
p0=[5 5]';	% Initial guess
% Old method:
% p=fminsearch(@dist2points, p0, [], points);
% New method:
p=fminsearch(@(x) dist2points(x, points), p0);
fprintf('p=%s\n', mat2str(p));

% ====== Create the objective function surface
x=linspace(-1, 4, 51);
y=linspace(-1, 5, 51);
zz=zeros(length(y), length(x));
for i=1:length(x)
	for j=1:length(y);
		zz(j,i)=dist2points([x(i), y(j)]', points);
	end
end
% ====== Mesh & contour plot
subplot(121); h=meshc(x, y, zz); axis image; box on
xlabel('X'); ylabel('Y'); zlabel('Z');
set(h(1), 'faceAlpha', 0.1);	% Make the mesh transparent
axisLimit=axis;
line(p(1), p(2), axisLimit(5), 'marker', '.', 'color', 'r');
for i=1:size(points,2)
	line(points(1,i), points(2,i), axisLimit(5)*[1 1 1 1], 'marker', '.', 'color', 'k', 'markerSize', 20);
end

% ====== Detailed contour plot
subplot(122); contour(x, y, zz, 100); axis image
xlabel('X'); ylabel('Y'); zlabel('Z');
line(p(1), p(2), 'marker', '.', 'color', 'r', 'markerSize', 20);
for i=1:size(points,2)
	line(points(1,i), points(2,i), 'marker', '.', 'color', 'k', 'markerSize', 20);
end