% Script to generate results
% Type:
%   Results
% to run this script.

close all;

% Image are blurred by 3 different PSF's and then additive noise is
% injected. Two types of noise is used: Gaussian and Uniform, both
% additive.

blur = {'motion'; 'gaussian'; 'average'};
noise = [0.001, 0.01, 0.05];
filename = {'images/wheel.png'; 'images/camera.png'; 'images/aerial.png'};

%_______________________________________________________________________________
% This section generates images using Gaussian Noise.

b = blur{1}
for x=1:3
    I = double(imread(filename{x}));
    I = I ./ 255;
    PSF = Fspecial('motion', 10, -45);
    for w=1:3
        N = noise(w) .* randn(size(I));
        D = imadd(imfilter(I, PSF, 'circular'), N);
        W = wienerfilter(D, PSF, I, N, 1);
        inv = inversefilter(D, PSF);
        mkdir(sprintf('%s_gaussian_noise/', filename{x}));
        D = (D + ones(size(D))*-min(min(D))) ./ (max(max(D + ones(size(D))*-min(min(D)))));
        W = (W + ones(size(W))*-min(min(W))) ./ (max(max(W + ones(size(W))*-min(min(W)))));     
        inv = (inv + ones(size(inv))*-min(min(inv))) ./ (max(max(inv + ones(size(inv))*-min(min(inv)))));
        imwrite(D, sprintf('%s_gaussian_noise/degraded_image_%s_%d.png', filename{x}, b, noise(w)*1000), 'png');
        imwrite(W, sprintf('%s_gaussian_noise/Wfilter_image_%s_%d.png', filename{x}, b, noise(w)*1000), 'png');
        imwrite(inv, sprintf('%s_gaussian_noise/Invfilter_image_%s_%d.png', filename{x}, b, noise(w)*1000), 'png');
    end
end    

b = blur{2}
for x=1:3
    I = double(imread(filename{x}));
    I = I ./ 255;
    PSF = Fspecial('gaussian', [11 11], 0.9);
    for w=1:3
        N = noise(w) .* randn(size(I));
        D = imadd(imfilter(I, PSF, 'circular'), N);
        W = wienerfilter(D, PSF, I, N, 1);
        inv = inversefilter(D, PSF);    
        mkdir(sprintf('%s_gaussian_noise/', filename{x}));
        inv = (inv + ones(size(inv))*-min(min(inv))) ./ (max(max(inv + ones(size(inv))*-min(min(inv)))));
        D = (D + ones(size(D))*-min(min(D))) ./ (max(max(D + ones(size(D))*-min(min(D)))));
        W = (W + ones(size(W))*-min(min(W))) ./ (max(max(W + ones(size(W))*-min(min(W)))));     
        imwrite(D, sprintf('%s_gaussian_noise/degraded_image_%s_%d.png', filename{x}, b, noise(w)*1000), 'png');
        imwrite(W, sprintf('%s_gaussian_noise/Wfilter_image_%s_%d.png', filename{x}, b, noise(w)*1000), 'png');
        imwrite(inv, sprintf('%s_gaussian_noise/Invfilter_image_%s_%d.png', filename{x}, b, noise(w)*1000), 'png');        
    end
end    

b = blur{3}
for x=1:3
    I = double(imread(filename{x}));
    I = I ./ 255;
    PSF = Fspecial('average', [7 7]);
    for w=1:3
        N = noise(w) .* randn(size(I));
        D = imadd(imfilter(I, PSF, 'circular'), N);
        W = wienerfilter(D, PSF, I, N, 1);
        inv = inversefilter(D, PSF);
        mkdir(sprintf('%s_gaussian_noise/', filename{x}));
        D = (D + ones(size(D))*-min(min(D))) ./ (max(max(D + ones(size(D))*-min(min(D)))));
        W = (W + ones(size(W))*-min(min(W))) ./ (max(max(W + ones(size(W))*-min(min(W)))));     
        inv = (inv + ones(size(inv))*-min(min(inv))) ./ (max(max(inv + ones(size(inv))*-min(min(inv)))));
        imwrite(D, sprintf('%s_gaussian_noise/degraded_image_%s_%d.png', filename{x}, b, noise(w)*1000), 'png');
        imwrite(W, sprintf('%s_gaussian_noise/Wfilter_image_%s_%d.png', filename{x}, b, noise(w)*1000), 'png');
        imwrite(inv, sprintf('%s_gaussian_noise/Invfilter_image_%s_%d.png', filename{x}, b, noise(w)*1000), 'png');        
    end
end    
%_______________________________________________________________________________
% This section generates images using Uniform Noise.

b = blur{1}
for x=1:3
    I = double(imread(filename{x}));
    I = I ./ 255;
    PSF = Fspecial('motion', 10, -45);
    for w=1:3
        N = noise(w) .* rand(size(I));
        D = imadd(imfilter(I, PSF, 'circular'), N);
        W = wienerfilter(D, PSF, I, N, 1);
        inv = inversefilter(D, PSF);
        mkdir(sprintf('%s_uniform_noise/', filename{x}));
        D = (D + ones(size(D))*-min(min(D))) ./ (max(max(D + ones(size(D))*-min(min(D)))));
        W = (W + ones(size(W))*-min(min(W))) ./ (max(max(W + ones(size(W))*-min(min(W)))));     
        inv = (inv + ones(size(inv))*-min(min(inv))) ./ (max(max(inv + ones(size(inv))*-min(min(inv)))));
        imwrite(D, sprintf('%s_uniform_noise/degraded_image_%s_%d.png', filename{x}, b, noise(w)*1000), 'png');
        imwrite(W, sprintf('%s_uniform_noise/Wfilter_image_%s_%d.png', filename{x}, b, noise(w)*1000), 'png');
        imwrite(inv, sprintf('%s_uniform_noise/Invfilter_image_%s_%d.png', filename{x}, b, noise(w)*1000), 'png');
    end
end    

b = blur{2}
for x=1:3
    I = double(imread(filename{x}));
    I = I ./ 255;
    PSF = Fspecial('gaussian', [11 11], 0.9);
    for w=1:3
        N = noise(w) .* rand(size(I));
        D = imadd(imfilter(I, PSF, 'circular'), N);
        W = wienerfilter(D, PSF, I, N, 1);
        inv = inversefilter(D, PSF);    
        mkdir(sprintf('%s_uniform_noise/', filename{x}));
        inv = (inv + ones(size(inv))*-min(min(inv))) ./ (max(max(inv + ones(size(inv))*-min(min(inv)))));
        D = (D + ones(size(D))*-min(min(D))) ./ (max(max(D + ones(size(D))*-min(min(D)))));
        W = (W + ones(size(W))*-min(min(W))) ./ (max(max(W + ones(size(W))*-min(min(W)))));     
        imwrite(D, sprintf('%s_uniform_noise/degraded_image_%s_%d.png', filename{x}, b, noise(w)*1000), 'png');
        imwrite(W, sprintf('%s_uniform_noise/Wfilter_image_%s_%d.png', filename{x}, b, noise(w)*1000), 'png');
        imwrite(inv, sprintf('%s_uniform_noise/Invfilter_image_%s_%d.png', filename{x}, b, noise(w)*1000), 'png');        
    end
end    

b = blur{3}
for x=1:3
    I = double(imread(filename{x}));
    I = I ./ 255;
    PSF = Fspecial('average', [7 7]);
    for w=1:3
        N = noise(w) .* rand(size(I));
        D = imadd(imfilter(I, PSF, 'circular'), N);
        W = wienerfilter(D, PSF, I, N, 1);
        inv = inversefilter(D, PSF);
        mkdir(sprintf('%s_uniform_noise/', filename{x}));
        D = (D + ones(size(D))*-min(min(D))) ./ (max(max(D + ones(size(D))*-min(min(D)))));
        W = (W + ones(size(W))*-min(min(W))) ./ (max(max(W + ones(size(W))*-min(min(W)))));     
        inv = (inv + ones(size(inv))*-min(min(inv))) ./ (max(max(inv + ones(size(inv))*-min(min(inv)))));
        imwrite(D, sprintf('%s_uniform_noise/degraded_image_%s_%d.png', filename{x}, b, noise(w)*1000), 'png');
        imwrite(W, sprintf('%s_uniform_noise/Wfilter_image_%s_%d.png', filename{x}, b, noise(w)*1000), 'png');
        imwrite(inv, sprintf('%s_uniform_noise/Invfilter_image_%s_%d.png', filename{x}, b, noise(w)*1000), 'png');        
    end
end
%_______________________________________________________________________________