idamnjanovic@45: function reconstructed=ImgDenoise_reconstruct(y, Problem, SparseDict)
ivan@128: %%  Image Denoising Problem reconstruction function
ivan@128: %   
ivan@128: %   This reconstruction function is using sparse representation y 
ivan@128: %   in dictionary Problem.A to reconstruct the patches of the denoised
ivan@128: %   image.
ivan@128: 
idamnjanovic@45: %
idamnjanovic@45: %   Centre for Digital Music, Queen Mary, University of London.
idamnjanovic@45: %   This file copyright 2009 Ivan Damnjanovic.
idamnjanovic@45: %
idamnjanovic@45: %   This program is free software; you can redistribute it and/or
idamnjanovic@45: %   modify it under the terms of the GNU General Public License as
idamnjanovic@45: %   published by the Free Software Foundation; either version 2 of the
idamnjanovic@45: %   License, or (at your option) any later version.  See the file
idamnjanovic@45: %   COPYING included with this distribution for more information.
ivan@128: %%
idamnjanovic@45: 
idamnjanovic@45: 
idamnjanovic@45: % stepsize %
idamnjanovic@45: if (isfield(Problem,'stepsize'))
idamnjanovic@45:   stepsize = Problem.stepsize;
idamnjanovic@45:   if (numel(stepsize)==1)
idamnjanovic@45:     stepsize = ones(1,2)*stepsize;
idamnjanovic@45:   end
idamnjanovic@45: else
idamnjanovic@45:   stepsize = ones(1,2);
idamnjanovic@45: end
idamnjanovic@45: if (any(stepsize<1))
idamnjanovic@45:   error('Invalid step size.');
idamnjanovic@45: end
idamnjanovic@45: 
idamnjanovic@45: % lambda %
idamnjanovic@45: if (isfield(Problem,'lambda'))
idamnjanovic@45:   lambda = Problem.lambda;
idamnjanovic@45: else
idamnjanovic@45:   lambda = Problem.maxval/(10*Problem.sigma);
idamnjanovic@45: end
idamnjanovic@45: if exist('SparseDict','var')&&(SparseDict==1)
idamnjanovic@45:     if issparse(Problem.A)
idamnjanovic@45:         A = Problem.A;
idamnjanovic@45:       else
idamnjanovic@45:         A = sparse(Problem.A);
idamnjanovic@45:       end
idamnjanovic@45:     cl_samp=add_dc(dictsep(Problem.basedict,A,y), Problem.b1dc,'columns');
idamnjanovic@45: else
idamnjanovic@45:     cl_samp=add_dc(Problem.A*y, Problem.b1dc,'columns');
idamnjanovic@45: end
idamnjanovic@45: %   combine the patches into reconstructed image
idamnjanovic@45: cl_im=col2imstep(cl_samp, size(Problem.Noisy), Problem.blocksize);
idamnjanovic@45: 
idamnjanovic@45: cnt = countcover(size(Problem.Noisy),Problem.blocksize,stepsize);
idamnjanovic@45: 
idamnjanovic@45: im = (cl_im+lambda*Problem.Noisy)./(cnt + lambda);
idamnjanovic@45: % y(y~=0)=1;
idamnjanovic@45: % numD=sum(y,2);
idamnjanovic@45: % nnzy=sum(y,1);
idamnjanovic@45: % figure(200);plot(sort(numD));
idamnjanovic@45: % figure(201);plot(sort(nnzy));
ivan@116: [v.RMSErn, v.RMSEcd, v.rn_im, v.cd_im]=SMALL_vmrse_type2(Problem.Original, Problem.Noisy, im);
idamnjanovic@45: %% output structure image+psnr %%
idamnjanovic@45: reconstructed.Image=im;
idamnjanovic@45: reconstructed.psnr = 20*log10(Problem.maxval * sqrt(numel(Problem.Original(:))) / norm(Problem.Original(:)-im(:)));
idamnjanovic@45: reconstructed.vmrse=v;
ivan@116: reconstructed.ssim=SMALL_ssim_index(Problem.Original, im);
idamnjanovic@45: end