%%%% Pietro Perona -- October 8 2001
%%%% 
%%%%
%%%% Change the parameters NOISE and xy_std to reproduce the experiments done in class
%%%%


RR = 32; CC = 32; %% Image size
NOISE = 0.05;   %% Standard dev. of Gaussian noise added to each image
blank_ima = zeros(RR,CC);

L = 4; % Object size
OBJECT = ones(L,L);
N = 100; %% Number of training objects
xy_std = 10; %% Variability in the position of the object
xy = round(RR/2+xy_std*(rand(2,N)-0.5));


X = zeros(RR*CC,N);

for i=1:N,
    x = xy(1,i) + round([1:L]-L/2);
    y = xy(2,i) + round([1:L]-L/2);
    ima = blank_ima;
    ima(x,y) = OBJECT;
    ima = ima + NOISE*randn(size(ima));
    X(:,i) = ima(:);
    if i<10,
        figure(1); imagesc(ima); colormap(gray); drawnow; pause(0.5);
    end;
end;


[U,S,V] = svd(X);


figure(3); semilogy(diag(S));
figure(2);
for k=1:12,
    subplot(3,4,k); imagesc(reshape(U(:,k),RR,CC)); colormap(gray);
end;