图像质量评价

最新推荐文章于 2024-08-19 10:41:42 发布

yayan01

最新推荐文章于 2024-08-19 10:41:42 发布

阅读量6.6k

点赞数 4

分类专栏：标题

标题专栏收录该内容

3 篇文章 0 订阅

订阅专栏

在图像处理算法研究中，很多时候需要有客观评价指标来对算法的性能进行评价。

比如，在图像复原、图像滤波算法研究中，需要采用客观评价指标来定量的来测试算法恢复出的图像相对于参考图像的好坏程度。

本文介绍文献中提到到三个比较好的客观评价指标——峰值性噪比PSNR、模糊系数K、质量因素Q，其定义分别是：

这三个指标的详细定义见参考文献[1]~[3]，下面给出这三个评价指标的MatLab实现。

[html]view plaincopy 
    
 %说明：本文件为计算两幅视频图象相对于高清晰图象的质量，其中：  
   
 %eyechart1.bmp为未处理前质量较差图象，核心区域的截图保存为area_eyechart1.bmp  
   
 %eyechart2.bmp为某种算法处理后质量较好图象，核心区域的截图保存为area_eyechart2.bmp  
   
 %eyechart3.bmp为高清晰参考图象，核心区域的截图保存为area_eyechart3.bmp  
   
 %程序流程为  
   
 %第一步：先分别从eyechart1.bmp、eyechart2.bmp、eyechart3.bmp中截取出核心区域，并分别保存为area_eyechart1.bmp、area_eyechart1.bmp、area_eyechart3.bmp  
   
 %第二步：以area_eyechart3.bmp为参考图象，计算area_eyechart1.bmp的PSNR、模糊系数KBlur、质量指数Q  
   
 %第三步：以area_eyechart3.bmp为参考图象，计算area_eyechart2.bmp的PSNR、模糊系数KBlur、质量指数Q  
   
 %程序可直接运行，运行结果为：  
   
 %1.保存并显示生成的截图文件area_eyechart1.bmp、area_eyechart1.bmp、area_eyechart3.bmp  
   
 %2、在控制台先显示第二步的计算结果，即area_eyechart1.bmp的三个质量指标，然后接着显示第三步的计算结果，即area_eyechart2.bmp的三个质量指标  
   
 %运行结果分析：area_eyechart2.bmp的PSNR和质量指数Q高，表明其质量较好，而area_eyechart1.bmp的模糊系数KBlur较大  
   
 %这是因为其有很多噪声被当着了边缘能量来计算，这也从一个方面说明模糊系数KBlur的应用具有局限性,但前者KBlur>1，后者<1，理论上只有模糊的情况下  
   
 %KBlur是<=1的，这也可根据KBlur与1的关系来判定图象收噪声污染的程度.  
   
 %******************从eyechart1.bmp,eyechart2.bmp两个文件中截取测试图象区域,可保证两图象截取的区域严格对准*****  
   
 a=imread('eyechart1.bmp','bmp');  
   
 b=a([203:396],[249:440]);  
   
 a=imread('eyechart2.bmp','bmp');  
   
 c=a([203:396],[249:440]);  
   
 %*******************从eyechart3.bmp中截取测试参考图象，截取部分需要进行缩放，使之与eyechart1.bmp,eyechart2.bmp截取部分大小匹配*******************************************************************  
   
 a=imread('eyechart3.bmp','bmp');  
   
 d=a([62:406],[60:395]);  
   
 e=imresize(d,[length(b(:,1)),length(b(1,:))], 'bicubic');%由于eyechart3.bmp和eyechart1.bmp,eyechart2.bmp比例不一样，这里要进行比例调整  
   
 imwrite(b,'area_eyechart1.bmp','bmp');  
   
 imwrite(c,'area_eyechart2.bmp','bmp');  
   
 imwrite(e,'area_eyechart3.bmp','bmp');  
   
 subplot(1,3,1);  
   
 imshow(e);  
   
 title('eyechart3.bmp截取部分，参考图象');  
   
 hold on;  
   
 subplot(1,3,2);  
   
 imshow(b);  
   
 title('eyechart1.bmp截取部分');  
   
 hold on;  
   
 subplot(1,3,3);  
   
 imshow(c);  
   
 title('eyechart2.bmp截取部分');  
   
 %*******************以下部分为计算截取图象area_eyechart1.bmp和area_eyechart1.bmp的PSNR、模糊系数、质量指数Q*  
   
 % 本文件功能为对计算污染图象相对于源图象的质量  
   
 clc;  
   
 clear;  
   
 PSNRenable=1;%PSNR计算使能，为0不计算，为1，计算  
   
 KBlurenable=1;%模糊系数KBlur计算使能，为0不计算，为1，计算  
   
 Qenable=1;%质量指数Q计算使能，为0不计算，为1，计算  
   
 for m=1:2  
   
 imsrcnamehead='area_eyechart3';%源图象文件名头  
   
 imsrcnameext='bmp';%源图象文件名扩展  
   
 if m==1 %以area_eyechart1.bmp为测试图象  
   
 imdstname=strcat('area_eyechart1','.',imsrcnameext);%污染图象文件名，可修改  
   
 elseif m==2%以area_eyechart2.bmp为测试图象  
   
   imdstname=strcat('area_eyechart2','.',imsrcnameext);%污染图象文件名，可修改  
   
 end  
   
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  
   
 iminfo=imfinfo(strcat(imsrcnamehead,'.',imsrcnameext));%源图象信息读取  
   
 imsrc=imread(strcat(imsrcnamehead,'.',imsrcnameext));%源图象读取  
   
 imdst=imread(imdstname,imsrcnameext);%污染图象读取  
   
 doubleimsrc=double(imsrc);%转换为浮点类型  
   
 doubleimdst=double(imdst);%转换为浮点类型  
   
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%源图象和污染图象读取  
   
 W=iminfo.Width;%图象度  
   
 H=iminfo.Height;%图象高  
   
 %///PSNR计算/  
   
 if PSNRenable==1  
   
 PSNR=0.0;%PSNR赋初值  
   
 for j=1:H  
   
   for i=1:W  
   
       PSNR=PSNR+double((doubleimsrc(j,i)-doubleimdst(j,i))*(doubleimsrc(j,i)-doubleimdst(j,i)));  
   
   end  
   
 end  
   
 PSNR=PSNR/W/H;  
   
 PSNR=10*log10(255*255/PSNR)  
   
 %PSNR计算完毕//  
   
 end  
   
 %///模糊系数KBlur计算/  
   
 if KBlurenable==1  
   
 Sin=0.0;%Sin赋初值  
   
 Sout=0.0;  
   
 for j=2:H-1  
   
   for i=2:W-1  
   
       t=doubleimsrc(j-1,i+1)+doubleimsrc(j+1,i-1)-doubleimsrc(j-1,i-1)-doubleimsrc(j+1,i+1);  
   
       if t<0 t=-t;  
   
       end  
   
       Sin=Sin+t;%源图象邻域边缘能量计算  
   
       t=doubleimdst(j-1,i+1)+doubleimdst(j+1,i-1)-doubleimdst(j-1,i-1)-doubleimdst(j+1,i+1);  
   
       if t<0 t=-t;  
   
       end  
   
       Sout=Sout+t;%污染图象邻域边缘能量计算  
   
   end  
   
 end  
   
 KBlur=Sout/Sin  
   
 end  
   
 %KBlur计算完毕  
   
 %///质量指数Q计算//  
   
 if Qenable==1  
   
 Q=0.0;%Q赋初值  
   
 Qnum=0;%图象以7X7块大小计算每块的Q，逐象素的移动块窗口，这里Qnum为块数量的计数  
   
 for j=4:H-3  
   
   for i=4:W-3  
   
       midsrc=0.0;  
   
       middst=0.0;  
   
       varsrc=0.0;  
   
       vardst=0.0;%源图象和污染图象块内的平均值和方差赋初值  
   
       varsrcdst=0.0;%源图象和污染图象块内的协方差赋初值  
   
       for n=-3:3  
   
           for m=-3:3  
   
               midsrc=midsrc+doubleimsrc(j+n,i+m);  
   
               middst=middst+doubleimdst(j+n,i+m);  
   
           end  
   
       end  
   
       midsrc=midsrc/49;  
   
       middst=middst/49;  
   
       %源图象和污染图象块内的平均值计算  
   
       for n=-3:3  
   
           for m=-3:3  
   
              varsrc=varsrc+(doubleimsrc(j+n,i+m)-midsrc)*(doubleimsrc(j+n,i+m)-midsrc);  
   
              vardst=vardst+(doubleimdst(j+n,i+m)-middst)*(doubleimdst(j+n,i+m)-middst);  
   
              varsrcdst=varsrcdst+(doubleimsrc(j+n,i+m)-midsrc)*(doubleimdst(j+n,i+m)-middst);  
   
           end  
   
       end  
   
       varsrc=varsrc/48;  
   
       vardst=vardst/48;  
   
       varsrcdst=varsrcdst/48;  
   
       if ((varsrc+vardst)*(midsrc*midsrc+middst*middst))~=0 %分母不为零的块才计算质量指数Q  
   
       Q=Q+4*varsrcdst*midsrc*middst/((varsrc+vardst)*(midsrc*midsrc+middst*middst));  
   
       %源图象和污染图象块内Q计算完毕  
   
       Qnum=Qnum+1;%块计数加1  
   
       end  
   
   end  
   
 end  
   
 Q=Q/Qnum  
   
 end  
   
 %质量指数Q计算完毕/  
   
 end  
   
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%图象质量计算完毕

参考文献

[1]袁飞，黄联芬，姚彦，视频质量客观评价技术研究，标准、检测与仪器，2007（3）：91-94

[2]黄文辉 ,陈仁雷 ,张家谋，数字视频图像质量客观测量方法的改进与实现，北京邮电大学学报，2005（4）：87-90

[3]Zhou Wang, Hamid R.Sheikh , Alan C. Bovik，Objective video quality assessment(Chapter 41 in The Handbook of Video Databases: Design and Applications)., CRC Press, 2003(1041-1078)

摘自：http://blog.csdn.net/droidpioneer/article/details/6612441

SSIM参数

一种衡量两幅图像相似度的新指标，其值越大越好，最大为1，

经常用到图像处理中，特别在图像去噪处理中在图像相似度评价上全面超越SNR（signal to noise ratio）和PSNR（peak signal to noise ratio）。

结构相似性理论认为，自然图像信号是高度结构化的，即像素间有很强的相关性，特别是空域中最接近的像素，这种相关性蕴含着视觉场景中物体结构的重要信息；HVS的主要功能是从视野中提取结构信息，可以用对结构信息的度量作为图像感知质量的近似。结构相似性理论是一种不同于以往模拟HVS低阶的组成结构的全新思想，与基于HVS特性的方法相比，最大的区别是自顶向下与自底向上的区别。这一新思想的关键是从对感知误差度量到对感知结构失真度量的转变。它没有试图通过累加与心理物理学简单认知模式有关的误差来估计图像质量，而是直接估计两个复杂结构信号的结构改变，从而在某种程度上绕开了自然图像内容复杂性及多通道去相关的问题。

作为结构相似性理论的实现，结构相似度指数从图像组成的角度将结构信息定义为独立于亮度、对比度的，反映场景中物体结构的属性，并将失真建模为亮度、对比度和结构三个不同因素的组合。用均值作为亮度的估计，标准差作为对比度的估计，协方差作为结构相似程度的度量。

（ http://baike.baidu.com/view/1559204.htm）

--------------------------------------------------------------------------------------------

结构相似性指标 (英文： structural similarity index， SSIM index)是一种用以衡量两张数位影像相似程度的指标。当两张影像其中一张为无失真影像，另一张为失真后的影像，二者的结构相似性可以看成是失真影像的影像品质衡量指标。相较于传统所使用的影像品质衡量指标，像是峰值信噪比 (英文： PSNR )，结构相似性在影像品质的衡量上更能符合人眼对影像品质的判断。

实际使用时，简化起见，一般会将参数设为 $\alpha=\beta=\gamma=1$ 及 $C_3=C_2/2$ ，得到：

$\text{SSIM}(\mathbf{x},\mathbf{y})=\frac{(2\mu_x\mu_y+C_1)(2\sigma_{xy}+C_2)}{(\mu_x^2+\mu_y^2+C_1)(\sigma_x^2+\sigma_y^2+C_2)}$ 。

在计算两张影像的结构相似性指标时，会开一个局部性的视窗，一般为 $N$ × $N$ 的小区块，计算出视窗内信号的结构相似性指标，每次以像素为单位移动视窗，直到整张影像每个位置的局部结构相似性指标都计算完毕。将全部的局部结构相似性指标平均起来即为两张影像的结构相似性指标。

摘自：http://blog.sina.com.cn/s/blog_70a336750101m508.html

ssim算法的matlab代码

---转载 http://bbs.blackbap.org/thread-2742-1-1.html
function [mssim, ssim_map] = ssim(img1, img2, K, window, L)

%
%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
%             depends on the window size and the downsampling factor.
%
�sic Usage:
%    Given 2 test images img1 and img2, whose dynamic range is 0-255
%
%    [mssim, ssim_map] = ssim(img1, img2);
%
�vanced Usage:
%    User defined parameters. For example
%
%    K = [0.05 0.05];
%    window = ones(8);
%    L = 100;
%    [mssim, ssim_map] = ssim(img1, img2, K, window, L);
%
%Visualize 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

img1 = double(img1);
img2 = double(img2);
% automatic downsampling
f = max(1,round(min(M,N)/256));
%downsampling by f
%use a simple low-pass filter
if(f>1)
   lpf = ones(f,f);
   lpf = lpf/sum(lpf(:));
   img1 = imfilter(img1,lpf,'symmetric','same');
   img2 = imfilter(img2,lpf,'symmetric','same');
   img1 = img1(1:f:end,1:f:end);
   img2 = img2(1:f:end,1:f:end);
end
C1 = (K(1)*L)^2;
C2 = (K(2)*L)^2;
window = window/sum(sum(window));
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);

摘自：http://blog.sina.com.cn/s/blog_6c41e2f301013ovs.html