view Problems/ImgDenoise_reconstruct.m @ 193:cc540df790f4 danieleb

Simple example that demonstrated dictionary learning... to be completed

function reconstructed=ImgDenoise_reconstruct(y, Problem, SparseDict)
%%  Image Denoising Problem reconstruction function
%   
%   This reconstruction function is using sparse representation y 
%   in dictionary Problem.A to reconstruct the patches of the denoised
%   image.

%
%   Centre for Digital Music, Queen Mary, University of London.
%   This file copyright 2009 Ivan Damnjanovic.
%
%   This program is free software; you can redistribute it and/or
%   modify it under the terms of the GNU General Public License as
%   published by the Free Software Foundation; either version 2 of the
%   License, or (at your option) any later version.  See the file
%   COPYING included with this distribution for more information.
%%


% stepsize %
if (isfield(Problem,'stepsize'))
  stepsize = Problem.stepsize;
  if (numel(stepsize)==1)
    stepsize = ones(1,2)*stepsize;
  end
else
  stepsize = ones(1,2);
end
if (any(stepsize<1))
  error('Invalid step size.');
end

% lambda %
if (isfield(Problem,'lambda'))
  lambda = Problem.lambda;
else
  lambda = Problem.maxval/(10*Problem.sigma);
end
if exist('SparseDict','var')&&(SparseDict==1)
    if issparse(Problem.A)
        A = Problem.A;
      else
        A = sparse(Problem.A);
      end
    cl_samp=add_dc(dictsep(Problem.basedict,A,y), Problem.b1dc,'columns');
else
    cl_samp=add_dc(Problem.A*y, Problem.b1dc,'columns');
end
%   combine the patches into reconstructed image
cl_im=col2imstep(cl_samp, size(Problem.Noisy), Problem.blocksize);

cnt = countcover(size(Problem.Noisy),Problem.blocksize,stepsize);

im = (cl_im+lambda*Problem.Noisy)./(cnt + lambda);
% y(y~=0)=1;
% numD=sum(y,2);
% nnzy=sum(y,1);
% figure(200);plot(sort(numD));
% figure(201);plot(sort(nnzy));
[v.RMSErn, v.RMSEcd, v.rn_im, v.cd_im]=SMALL_vmrse_type2(Problem.Original, Problem.Noisy, im);
%% output structure image+psnr %%
reconstructed.Image=im;
reconstructed.psnr = 20*log10(Problem.maxval * sqrt(numel(Problem.Original(:))) / norm(Problem.Original(:)-im(:)));
reconstructed.vmrse=v;
reconstructed.ssim=SMALL_ssim_index(Problem.Original, im);
end