# HG changeset patch # User idamnjanovic # Date 1300116798 0 # Node ID 84b7c44e54c4caac4bc8ef8c248fe7ffd49912d4 # Parent d80c103d9876954b79de52928e2844267fa2eefb diff -r d80c103d9876 -r 84b7c44e54c4 util/ssim_index.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/util/ssim_index.m Mon Mar 14 15:33:18 2011 +0000 @@ -0,0 +1,194 @@ +function [mssim, ssim_map] = ssim_index(img1, img2, K, window, L) + +%======================================================================== +%SSIM Index, Version 1.0 +%Copyright(c) 2003 Zhou Wang +%All Rights Reserved. +% +%The author is with Howard Hughes Medical Institute, and Laboratory +%for Computational Vision at Center for Neural Science and Courant +%Institute of Mathematical Sciences, New York University. +% +%---------------------------------------------------------------------- +%Permission to use, copy, or modify this software and its documentation +%for educational and research purposes only and without fee is hereby +%granted, provided that this copyright notice and the original authors' +%names appear on all copies and supporting documentation. This program +%shall not be used, rewritten, or adapted as the basis of a commercial +%software or hardware product without first obtaining permission of the +%authors. The authors make no representations about the suitability of +%this software for any purpose. It is provided "as is" without express +%or implied warranty. +%---------------------------------------------------------------------- +% +%This is an implementation of the algorithm for calculating the +%Structural SIMilarity (SSIM) index between two images. Please refer +%to the following paper: +% +%Z. Wang, A. C. Bovik, H. R. Sheikh, and E. P. Simoncelli, "Image +%quality assessment: From error measurement to structural similarity" +%IEEE Transactios on Image Processing, vol. 13, no. 1, Jan. 2004. +% +%Kindly report any suggestions or corrections to zhouwang@ieee.org +% +%---------------------------------------------------------------------- +% +%Input : (1) img1: the first image being compared +% (2) img2: the second image being compared +% (3) K: constants in the SSIM index formula (see the above +% reference). defualt value: K = [0.01 0.03] +% (4) window: local window for statistics (see the above +% reference). default widnow is Gaussian given by +% window = fspecial('gaussian', 11, 1.5); +% (5) L: dynamic range of the images. default: L = 255 +% +%Output: (1) mssim: the mean SSIM index value between 2 images. +% If one of the images being compared is regarded as +% perfect quality, then mssim can be considered as the +% quality measure of the other image. +% If img1 = img2, then mssim = 1. +% (2) ssim_map: the SSIM index map of the test image. The map +% has a smaller size than the input images. The actual size: +% size(img1) - size(window) + 1. +% +%Default Usage: +% Given 2 test images img1 and img2, whose dynamic range is 0-255 +% +% [mssim ssim_map] = ssim_index(img1, img2); +% +%Advanced Usage: +% User defined parameters. For example +% +% K = [0.05 0.05]; +% window = ones(8); +% L = 100; +% [mssim ssim_map] = ssim_index(img1, img2, K, window, L); +% +%See the results: +% +% mssim %Gives the mssim value +% imshow(max(0, ssim_map).^4) %Shows the SSIM index map +% +%======================================================================== + + +if (nargin < 2 || nargin > 5) + ssim_index = -Inf; + ssim_map = -Inf; + return; +end + +if (size(img1) ~= size(img2)) + ssim_index = -Inf; + ssim_map = -Inf; + return; +end + +[M N] = size(img1); + +if (nargin == 2) + if ((M < 11) || (N < 11)) + ssim_index = -Inf; + ssim_map = -Inf; + return + end + window = fspecial('gaussian', 11, 1.5); % + K(1) = 0.01; % default settings + K(2) = 0.03; % + L = 255; % +end + +if (nargin == 3) + if ((M < 11) || (N < 11)) + ssim_index = -Inf; + ssim_map = -Inf; + return + end + window = fspecial('gaussian', 11, 1.5); + L = 255; + if (length(K) == 2) + if (K(1) < 0 || K(2) < 0) + ssim_index = -Inf; + ssim_map = -Inf; + return; + end + else + ssim_index = -Inf; + ssim_map = -Inf; + return; + end +end + +if (nargin == 4) + [H W] = size(window); + if ((H*W) < 4 || (H > M) || (W > N)) + ssim_index = -Inf; + ssim_map = -Inf; + return + end + L = 255; + if (length(K) == 2) + if (K(1) < 0 || K(2) < 0) + ssim_index = -Inf; + ssim_map = -Inf; + return; + end + else + ssim_index = -Inf; + ssim_map = -Inf; + return; + end +end + +if (nargin == 5) + [H W] = size(window); + if ((H*W) < 4 || (H > M) || (W > N)) + ssim_index = -Inf; + ssim_map = -Inf; + return + end + if (length(K) == 2) + if (K(1) < 0 || K(2) < 0) + ssim_index = -Inf; + ssim_map = -Inf; + return; + end + else + ssim_index = -Inf; + ssim_map = -Inf; + return; + end +end + +C1 = (K(1)*L)^2; +C2 = (K(2)*L)^2; +window = window/sum(sum(window)); +img1 = double(img1); +img2 = double(img2); + +mu1 = filter2(window, img1, 'valid'); +mu2 = filter2(window, img2, 'valid'); +mu1_sq = mu1.*mu1; +mu2_sq = mu2.*mu2; +mu1_mu2 = mu1.*mu2; +sigma1_sq = filter2(window, img1.*img1, 'valid') - mu1_sq; +sigma2_sq = filter2(window, img2.*img2, 'valid') - mu2_sq; +sigma12 = filter2(window, img1.*img2, 'valid') - mu1_mu2; + +if (C1 > 0 & C2 > 0) + ssim_map = ((2*mu1_mu2 + C1).*(2*sigma12 + C2))./((mu1_sq + mu2_sq + C1).*(sigma1_sq + sigma2_sq + C2)); +else + numerator1 = 2*mu1_mu2 + C1; + numerator2 = 2*sigma12 + C2; + denominator1 = mu1_sq + mu2_sq + C1; + denominator2 = sigma1_sq + sigma2_sq + C2; + ssim_map = ones(size(mu1)); + index = (denominator1.*denominator2 > 0); + ssim_map(index) = (numerator1(index).*numerator2(index))./(denominator1(index).*denominator2(index)); + index = (denominator1 ~= 0) & (denominator2 == 0); + ssim_map(index) = numerator1(index)./denominator1(index); +end + +mssim = mean2(ssim_map); + +return \ No newline at end of file diff -r d80c103d9876 -r 84b7c44e54c4 util/vmrse_type2.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/util/vmrse_type2.m Mon Mar 14 15:33:18 2011 +0000 @@ -0,0 +1,43 @@ +function [RMSErn, RMSEcd, rn_im, cd_im] = vmrse_type2(orig, corr, recon) + +%%% Implementation of VectorRMSE type2 +% +% Centre for Digital Music, Queen Mary, University of London. +% This file copyright 2011 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. +% +% Input: +% - Original image +% - Corrupted image +% - Reconstructed Image +% +% Output: +% - RMSErn - RMSE from residual noise (noise not completely removed) +% - RMSEcd - RMSE from collateral distortion - excessive filtering +% - rn_im - image of residual noise +% - cd_im - image of collateral distortion +% +% F. Russo, "New Method for Performance Evaluation of Grayscale Image +% Denoising filters", IEEE Signal Processing Letters, vol. 17, no. 5, +% pp.417-420, May 2010 +%% + + recon_int = round(recon); + + RN1 = ((origrecon_int)&(recon_int>=corr)); + CD1 = ((origcorr)); + CD2 = ((orig>recon_int)&(recon_int