function [G, e] = gersh(A, noplot)
%GERSH Gershgorin disks.
% GERSH(A) draws the Gershgorin disks for the square matrix A.
% The eigenvalues are plotted as crosses `x'.
% Alternative usage: [G, E] = GERSH(A, 1) suppresses the plot
% and returns the data in G, with A's eigenvalues in E.
%
% Try GERSH(GALLERY('LESP',N)) and GERSH(GALLERY('SMOKE',N)).
if diff(size(A)), error('Matrix must be square.'), end
n = length(A);
m = 40;
G = zeros(m,n);
d = diag(A);
r = sum( abs( A-diag(d) )' )';
e = eig(A);
radvec = exp(i * linspace(0,2*pi,m)');
for j=1:n
G(:,j) = d(j)*ones(m,1) + r(j)*radvec;
end
if nargin < 2
ax = cpltaxes(G(:));
for j=1:n
plot(real(G(:,j)), imag(G(:,j)),'-') % Plot the disks.
hold on
end
plot(real(e), imag(e), 'x') % Plot the eigenvalues too.
axis(ax)
axis('square')
hold off
end