Matlab常见的图像质量评价指标

(一)有参考图像评价指标
1.信噪比SNR

function SNR=imSNR(In,I)
% 计算信号噪声比函数
% I :original signal
% In:noisy signal(ie. Original signal + noise signal)
% snr=10*log10(sigma2(I2)/sigma2(I2-I1))
[row,col,nchannel]=size(I);
In=double(In);
I=double(I);
msnr=0;
if nchannel==1%gray image
    Ps=sum(sum((I-mean(mean(I))).^2));%signal power
    Pn=sum(sum((I-In).^2));%noise power
    msnr=10*log10(Ps/Pn);
elseif nchannel==3%color image
    for i=1:3
        Ps=sum(sum((I(:,:,i)-mean(mean(I(:,:,i)))).^2));%signal power
        Pn=sum(sum((I(:,:,i)-In(:,:,i)).^2));%noise power
        msnr=msnr+10*log10(Ps/Pn);
    end
    msnr=msnr/3;
end
imshow(uint8(In));
title(sprintf(' SNR :%.3f dB ', msnr));
return

2.峰值信噪比PSNR

function[PSNR]=imPSNR(J,I)
%imPSNR Summary of this function goes here
%I is a image with high quality
%J is a image with noise
%the function will return the PSNR of the noise image
if any(size(J)~=size(I))
    disp('Please check the input image have the same size');
end
J=double(J);
img=double(I);
[h,w]=size(I);
B=8;%编码一个像素用多少二进制位   
MAX=2^B-1;%图像有多少灰度级 
MSE=sum(sum((J-img).^2))/(h*w);%均方差   
PSNR=20*log10(MAX/sqrt(MSE));%峰值信噪比   
imshow(uint8(J));
title(sprintf(' PSNR :%.3f dB ', PSNR));

3.结构相似性SSIM

function [mssim, ssim_map,siga_sq,sigb_sq] = SSIM(ima, imb)  
%http://blog.csdn.net/cool_ben/article/details/46597967
% ========================================================================  
%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为灰度级数，L=255  
%-------------------------------------------------------------------  
%     ima - 比较图像A  
%     imb - 比较图像B  
%  
% ssim_map - 各加窗后得到的SSIM（A,B|w）组成的映射矩阵  
%    mssim - 对加窗得到的SSIM（A,B|w）求平均，即最终的SSIM（A,B）  
%  siga_sq - 图像A各窗口内灰度值的方差  
%  sigb_sq - 图像B各窗口内灰度值的方差  
%-------------------------------------------------------------------  
%  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');%对窗口内并没有进行平均处理，而是与高斯卷积，  
ub=filter2(w, imb, 'valid'); % 类似加权平均  
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);  
return 

4.相似性corr2

BPoutputSg2=imread('E:\matlab code\BPNN\BPNN_edge\moreFrameImg\allPic\BPimg\BPoutputSg2.bmp');
BPoutputSg10=imread('E:\matlab code\BPNN\BPNN_edge\moreFrameImg\allPic\BPimg\BPoutputSg10.bmp');
trueImg=imread('E:\matlab code\BPNN\BPNN_edge\moreFrameImg\allPic\BPimg\trueImg.bmp');
BPoutputSg2=BPoutputSg2(:,:,1);
BPoutputSg10=BPoutputSg10(:,:,1);
trueImg=trueImg(:,:,1);

r10=corr2(BPoutputSg10,trueImg);
imshow(BPoutputSg10);
title(sprintf('coefficient of association :%.4f', r10));
figure;
r2=corr2(BPoutputSg2,trueImg);
imshow(BPoutputSg2);
title(sprintf('coefficient of association :%.4f', r2));

(二)无参考图像评价指标
1.平均梯度

function AVEGRAD=avegrad(img)
%%%%%%%%%%%%%%%%%%%%%%%%%
%%%% this function is used to calculate
%%%% the average gradient of an image.
%%%% editor: denghaibo Email:denghblzu07@gmail.com
%%%% date : 2009-4-15
%%%% 平均梯度可敏感地反映图像对微小细节反差表达的能力，可用来评价图像的模糊程度
%%%% 在图像中，某一方向的灰度级变化率大，它的梯度也就大。因此，可以用平均梯度值
%%%% 来衡量图像的清晰度，还同时反映出图像中微小细节反差和纹理变换特征。

img=double(img);
[M,N]=size(img);
gradval=zeros(M,N); %%% save the gradient value of single pixel
diffX=zeros(M,N);    %%% save the differential value of X orient
diffY=zeros(M,N);    %%% save the differential value of Y orient

tempX=zeros(M,N);
tempY=zeros(M,N);
tempX(1:M,1:(N-1))=img(1:M,2:N);
tempY(1:(M-1),1:N)=img(2:M,1:N);

diffX=img-tempX;
diffY=img-tempY;
diffX(1:M,N)=0;       %%% the boundery set to 0
diffY(M,1:N)=0;
diffX=diffX.*diffX;
diffY=diffY.*diffY;
AVEGRAD=sum(sum(diffX+diffY));
AVEGRAD=sqrt(AVEGRAD/2);
AVEGRAD=AVEGRAD/((M-1)*(N-1));
disp(AVEGRAD);

end 

2.Brenner梯度

 [row,col]=size(img_after);
Bre=0;
for j=1:col
    for i=1:row-2
        Bre=Bre+(img_after(i+2,j)-img_after(i,j)).^2;
    end
end

3.图像熵Entropy
可以直接利用entropy函数实现

4.能量梯度函数EOG

 [row,col]=size(img_after);
EOG=0;
for i=1:row-1
    for j=1:col-1
        EOG=EOG+(img_after(i,j+1)-img_after(i,j)).^2+(img_after(i+1,j)-img_after(i,j)).^2;
    end
end

5.空间频率

[row,col]=size(img_after);
RF1=0;
RF1_sin=0;
for i=1:row
    for j=2:col
        RF1=RF1+(img_after(i,j)-img_after(i,j-1)).^2;
    end
end
RF=sqrt(1/(row*col)*RF1);

CF1=0;
for i=2:row
    for j=1:col
        CF1= CF1+(img_after(i,j)-img_after(i-1,j)).^2;
    end
end
CF=sqrt(1/(row*col)*CF1);
SF=sqrt(RF.^2+CF.^2);

6.方差或标准差

 [row,col]=size(img_after);
sum=0;
for i=1:row
    for j=1:col
        sum=sum+img_after(i,j);
    end
end
ave=(1/(row*col))*sum;
Var1=0;
for i=1:row
    for j=1:col
        Var1=Var1+(img_after(i,j)-ave).^2;
    end
end
Var=(1/(row*col))*Var1;
SD=sqrt(Var);