% Read in the image
[X,map]=imread('Fractal01.png');
image(X)
colormap(map)
axis square

%Force it to be purely symmetric
A=zeros(600,600);
for j=1:600
    A(j,:)=X(j,:);
    A(:,j)=X(j,:)';
end
image(A)
axis square

%Find the eigenvalues (in matrix D) and eigenvecs (in V).
[V,D]=eig(A);

%Plot the spectral decomposition, ordered by the magnitude of the
%eigenvalues:
t=diag(D);
[val,idx]=sort(abs(t),'descend');

for j=1:100
    H=V(:,idx(1:j))*diag(t(idx(1:j)))*V(:,idx(1:j))';
    image(H)
    Err(j)=norm(H-A,'fro');
    pause(0.1)
end