圖像相似度評價(jià)指標(biāo)
在圖像處理中我們經(jīng)常遇到需要評價(jià)兩張圖像是否相似做鹰，給出其相似度的指標(biāo)立莉，這里總結(jié)了三種評判指標(biāo)均方誤差MSE, 結(jié)構(gòu)相似性SSIM, 以及峰值信噪比PSNR, 分三個(gè)小結(jié)介紹其原理以及對應(yīng)的matlab以及tensorflow版本的算法實(shí)現(xiàn)昨悼。
均方誤差MSE
即m×n單色圖像 I 和 K（原圖像與處理圖像）之間均方誤差捣辆，定義為：
MSE formula
Matlab實(shí)現(xiàn)
function MSE = mse(I, K)
[M,N,D] = size(I);
Diff = double(I)-double(K);
MSE = sum(Diff(:).^2)/numel(I);
end
Tensorflow實(shí)現(xiàn)
# import tensorflow as tf
def MSE(I, K):
x, y = tf.cast(I, tf.float32), tf.cast(K, tf.float32)
mse = tf.losses.mean_squared_error(labels=y, predictions=x)
return mse
結(jié)構(gòu)相似性SSIM
- 結(jié)構(gòu)相似性：
自然圖像具有極高的結(jié)構(gòu)性换况，表現(xiàn)在圖像的像素間存在著很強(qiáng)的相關(guān)性欠气，尤其是在空間相似的情況下厅各。這些相關(guān)性在視覺場景中攜帶著關(guān)于物體結(jié)構(gòu)的重要信息。我們假設(shè)人類視覺系統(tǒng)（HSV）主要從可視區(qū)域內(nèi)獲取結(jié)構(gòu)信息预柒。所以通過探測結(jié)構(gòu)信息是否改變來感知圖像失真的近似信息队塘。
大多數(shù)的基于誤差敏感度（error sensitivity）的質(zhì)量評估方法(如MSE,PSNR)使用線性變換來分解圖像信號，這不會涉及到相關(guān)性宜鸯。我們要討論的SSIM就是要找到更加直接的方法來比較失真圖像和參考圖像的結(jié)構(gòu)憔古。- SSIM指數(shù)
物體表面的亮度信息與照度和反射系數(shù)有關(guān)，且場景中的物體的結(jié)構(gòu)與照度是獨(dú)立的顾翼，反射系數(shù)與物體有關(guān)投放。我們可以通過分離照度對物體的影響來探索一張圖像中的結(jié)構(gòu)信息。這里适贸，把與物體結(jié)構(gòu)相關(guān)的亮度和對比度作為圖像中結(jié)構(gòu)信息的定義灸芳。因?yàn)橐粋€(gè)場景中的亮度和對比度總是在變化的，所以我們可以通過分別對局部的處理來得到更精確的結(jié)果拜姿。
SSIM的算法流程圖原理圖如下所示：
SSIM測量系統(tǒng)
SSIM的求解公式如下：
SSIM formula
其中u_x是x的平均值烙样，u_y是y的平均值，σ_x是x的方差蕊肥，σ_y是y的方差谒获，σ_{xy}是x和y的協(xié)方差。c_1=(k_1*L)^2壁却，c_2=(k_2*L)^2是用來維持穩(wěn)定的常數(shù)批狱。L是像素值的動態(tài)范圍。k_1=0.01,k_2=0.03展东。
結(jié)構(gòu)相似性的范圍為-1到+1（即SSIM∈(-1, 0]）赔硫。當(dāng)兩張圖像一模一樣時(shí)，SSIM的值等于1盐肃。
Matlab實(shí)現(xiàn)
function [mssim, ssim_map,siga_sq,sigb_sq] = SSIM(ima, imb)
% ========================================================================
%ssim的算法主要參考如下論文：
%Z. Wang, A. C. Bovik, H. R. Sheikh, and E. P. Simoncelli, "Image
% quality assessment: From error visibility to structural similarity,"
% IEEE Transactios on Image Processing, vol. 13, no. 4, pp. 600-612,
% Apr. 2004.
% 首先對圖像加窗處理爪膊，w=fspecial('gaussian', 11, 1.5);
% (2*ua*ub+C1)*(2*sigmaa*sigmab+C2)
% SSIM(A,B)=——————————————————————————————————————————————————
% (ua*ua+ub*ub+C1)(sigmaa*sigmaa+sigmab*sigmab+C2)
% C1=（K1*L）;
% C2=(K2*L); K1=0.01，K2=0.03
% L為灰度級數(shù)砸王，L=255
%-------------------------------------------------------------------
% ima - 比較圖像A
% imb - 比較圖像B
%
% ssim_map - 各加窗后得到的SSIM（A,B|w）組成的映射矩陣
% mssim - 對加窗得到的SSIM（A,B|w）求平均推盛，即最終的SSIM（A,B）
% siga_sq - 圖像A各窗口內(nèi)灰度值的方差
% sigb_sq - 圖像B各窗口內(nèi)灰度值的方差
%-------------------------------------------------------------------
% Cool_ben
%=======================================================================
w = fspecial('gaussian', 11, 1.5); %window 加窗
K(1) = 0.01;
K(2) = 0.03;
L = 255;
ima = double(ima);
imb = double(imb);
C1 = (K(1)*L)^2;
C2 = (K(2)*L)^2;
w = w/sum(sum(w));
ua = filter2(w, ima, 'valid');%對窗口內(nèi)并沒有進(jìn)行平均處理，而是與高斯卷積谦铃，
ub = filter2(w, imb, 'valid'); % 類似加權(quán)平均
ua_sq = ua.*ua;
ub_sq = ub.*ub;
ua_ub = ua.*ub;
siga_sq = filter2(w, ima.*ima, 'valid') - ua_sq;
sigb_sq = filter2(w, imb.*imb, 'valid') - ub_sq;
sigab = filter2(w, ima.*imb, 'valid') - ua_ub;
ssim_map = ((2*ua_ub + C1).*(2*sigab + C2))./((ua_sq + ub_sq + C1).*(siga_sq + sigb_sq + C2));
mssim = mean2(ssim_map);
end
Tensorflow實(shí)現(xiàn)
# import tensorflow as tf
def _tf_fspecial_gauss(size, sigma):
"""Function to mimic the 'fspecial' gaussian MATLAB function"""
x_data, y_data = np.mgrid[-size//2 + 1:size//2 + 1, -size//2 + 1:size//2 + 1]
x_data = np.expand_dims(x_data, axis=-1)
x_data = np.expand_dims(x_data, axis=-1)
y_data = np.expand_dims(y_data, axis=-1)
y_data = np.expand_dims(y_data, axis=-1)
x = tf.constant(x_data, dtype=tf.float32)
y = tf.constant(y_data, dtype=tf.float32)
g = tf.exp(-((x**2 + y**2)/(2.0*sigma**2)))
return g / tf.reduce_sum(g)
def tf_ssim(img1, img2, cs_map=False, mean_metric=True, size=11, sigma=1.5):
window = _tf_fspecial_gauss(size, sigma) # window shape [size, size]
K1 = 0.01
K2 = 0.03
L = 1 # depth of image (255 in case the image has a different scale)
C1 = (K1*L)**2
C2 = (K2*L)**2
mu1 = tf.nn.conv2d(img1, window, strides=[1, 1, 1, 1], padding='VALID')
mu2 = tf.nn.conv2d(img2, window, strides=[1, 1, 1, 1], padding='VALID')
mu1_sq = mu1*mu1
mu2_sq = mu2*mu2
mu1_mu2 = mu1*mu2
sigma1_sq = tf.nn.conv2d(img1*img1, window, strides=[1, 1, 1, 1], padding='VALID') - mu1_sq
sigma2_sq = tf.nn.conv2d(img2*img2, window, strides=[1, 1, 1, 1], padding='VALID') - mu2_sq
sigma12 = tf.nn.conv2d(img1*img2, window, strides=[1, 1, 1, 1], padding='VALID') - mu1_mu2
if cs_map:
value = (((2*mu1_mu2 + C1)*(2*sigma12 + C2))/((mu1_sq + mu2_sq + C1) *
(sigma1_sq + sigma2_sq + C2)),
(2.0*sigma12 + C2)/(sigma1_sq + sigma2_sq + C2))
else:
value = ((2*mu1_mu2 + C1)*(2*sigma12 + C2))/((mu1_sq + mu2_sq + C1) *
(sigma1_sq + sigma2_sq + C2))
if mean_metric:
value = tf.reduce_mean(value)
return value
def tf_ms_ssim(img1, img2, mean_metric=True, level=5):
weight = tf.constant([0.0448, 0.2856, 0.3001, 0.2363, 0.1333], dtype=tf.float32)
mssim = []
mcs = []
for l in range(level):
ssim_map, cs_map = tf_ssim(img1, img2, cs_map=True, mean_metric=False)
mssim.append(tf.reduce_mean(ssim_map))
mcs.append(tf.reduce_mean(cs_map))
filtered_im1 = tf.nn.avg_pool(img1, [1, 2, 2, 1], [1, 2, 2, 1], padding='SAME')
filtered_im2 = tf.nn.avg_pool(img2, [1, 2, 2, 1], [1, 2, 2, 1], padding='SAME')
img1 = filtered_im1
img2 = filtered_im2
# list to tensor of dim D+1
mssim = tf.stack(mssim, axis=0)
mcs = tf.stack(mcs, axis=0)
value = (tf.reduce_prod(mcs[0:level-1]**weight[0:level-1])*(mssim[level-1]**weight[level-1]))
if mean_metric:
value = tf.reduce_mean(value)
return value
或者耘成，emm，Tensorflow r1.8考慮使用下面函數(shù)：
tf.image.ssim(
img1,
img2,
max_val
)
具體使用方式參考：傳送門
峰值信噪比PSNR
PSNR本質(zhì)上與MSE相同，是MSE的對數(shù)表示凿跳。
峰值信噪比PSNR衡量圖像失真或是噪聲水平的客觀標(biāo)準(zhǔn)件豌。2個(gè)圖像之間PSNR值越大，則越相似控嗜。普遍基準(zhǔn)為30dB，30dB以下的圖像劣化較為明顯骡显。定義為：
PSNR formula
Matlab實(shí)現(xiàn)
function [PSNR, MSE]=psnr(I,K)
[M,N,D] = size(I);
Diff = double(I)-double(K);
MSE = sum(Diff(:).^2)/numel(I);
PSNR=10*log10(255^2/MSE);
end
Tensorflow實(shí)現(xiàn)
# import tensorflow as tf
def PSNR(I, K):
x, y = tf.cast(I, tf.float32), tf.cast(K, tf.float32)
mse = tf.losses.mean_squared_error(labels=y, predictions=x)
psnr = 10*tf.log(255**2/mse)/tf.log(10)
return psnr
