% wiener2_template.m adaptive filtering template

rng(0)
nx = 62; ny = 50;

p = 0.004;
hh = ones(9,7) / sqrt(p*(1-p));
ff = conv2(double(rand(nx,ny) < p), hh, 'same') - p; % true image
gg = ff + 5 * randn(size(ff));	% add gaussian noise

f.w2 = wiener2(gg, [5 5], 25); % conventional adapative wiener
f.npls = npls_sps(gg);

% REPLACE THE FOLLOWING LINE WITH YOUR OWN IMPROVED METHOD!
f.new = wiener2(gg, [5 5], 40); % just a place-holder for script

clf, pl = @(i) subplot(330+i);
ilim = [0 55]; elim = [-10 10]; colormap(1-gray(256))
xy = @() set(gca, 'xtick', [1 nx], 'ytick', [1 ny]);

pl(1), imagesc(1:nx,1:ny,ff', ilim), axis xy, axis image
xlabel n, ylabel m, title 'True: f[n,m]', colorbar horiz, xy()

pl(2), imagesc(1:nx,1:ny,gg', ilim), axis xy, axis image
xlabel n, ylabel m, colorbar horiz, xy()
title(sprintf('Noisy g[n,m], RMSE=%.2f', sqrt(mean2((gg-ff).^2))))

pl(4), imagesc(1:nx, 1:ny, f.w2', ilim), axis xy, axis image
xlabel n, ylabel m, title 'Wiener2 output', colorbar horiz, xy()

pl(5), imagesc(1:nx, 1:ny, f.npls', ilim), axis xy, axis image
xlabel n, ylabel m, title 'NPLS output', colorbar horiz, xy()

pl(6), imagesc(1:nx, 1:ny, f.new', ilim), axis xy, axis image
xlabel n, ylabel m, title 'My Method!', colorbar horiz, xy()

pl(7), imagesc(1:nx, 1:ny, (f.w2-ff)', elim), axis xy, axis image
xlabel n, ylabel m, colorbar horiz, xy()
title(sprintf('Wiener2 error, RMSE=%.2f', sqrt(mean2((f.w2-ff).^2))))

pl(8), imagesc(1:nx, 1:ny, (f.npls-ff)', elim), axis xy, axis image
xlabel n, ylabel m, colorbar horiz, xy()
title(sprintf('NPLS error, RMSE=%.2f', sqrt(mean2((f.npls-ff).^2))))

pl(9), imagesc(1:nx, 1:ny, (f.new-ff)', elim), axis xy, axis image
xlabel n, ylabel m, colorbar horiz, xy()
title(sprintf('My Error, RMSE=%.2f', sqrt(mean2((f.new-ff).^2))))

%orient tall, print('wiener2p', '-deps')
