数字图像处理，叠加平均去噪与小波阈值去噪的对比

叠加平均去噪：反复对同一张干净的图像添加随机噪声（椒盐，乘性，高斯等噪声，此处噪声不叠加，具体见代码），进行200张噪声图像的平均。

小波阈值去噪：利用改进了阈值的小波阈值去噪对图像进行去噪处理，其中半阈值是改进了的阈值函数。

图像信息：经典的lena图，大小256*256，256位灰度图

噪声浓度均为0.05

1，添加'speckle'噪声时，两种算法的去砸效果对比

利用客观参数psnr评估可知：

pnsr_averge =40.9679
pnsr_noise =18.9210
pnsr_hard =23.9319
pnsr_soft =24.4425
pnsr_half =25.0121

结论1，清晰灰度图像下的乘性噪声污染，不管是主观还是客观参数，都说明叠加平均的去噪能力远优于小波去噪，叠加平均效果非常好。

2，添加'salt & pepper'噪声时，两种算法的去砸效果对比

pnsr_averge =38.0104
pnsr_noise =18.3476
pnsr_hard =20.2537
pnsr_soft =24.0672
pnsr_half =22.5696

结论2，清晰灰度图像下的椒盐噪声污染，不管是主观还是客观参数，都说明叠加平均的去噪能力远优于小波去噪，叠加平均效果非常好。

3，添加 'gaussian' 噪声时，两种算法的去砸效果对比

pnsr_averge =25.9784
pnsr_noise =19.0491
pnsr_hard =22.6575
pnsr_soft =22.1231
pnsr_half =23.0135

结论3：结论还是一致的，叠加平均效果非常好。

4，添加浓度（0.1）更高更复杂的混合（三种噪声的混合）噪声，两种算法的去砸效果对比

pnsr_averge =20.7570
pnsr_noise =12.5467
pnsr_hard =18.3766
pnsr_soft =18.7248
pnsr_half =18.4966

最终结论：

本文简单地以小波阈值去噪为辅助参考效果，并采用客观评价指标得出叠加平均去噪法对于随机噪声具有非常强的去噪能力。该方法适合对同一景物的多次采集图像相加后取平均值来消除噪声。在实际运算过程中不需要图像的统计特性,可以克服线性滤波器所带来的图像细节模糊,而且对随机的滤波脉冲干扰和散斑图像噪声最有效。设想，如果同一景物的多次采集图像由于不可抗拒因素存在微小的位移，那么此时就需要用配准来校准了。

参考代码：

clc;
clear;
close all;

src_img=imread('lena.bmp');
%I=rgb2gray(I);
figure(1);
subplot(231)
imshow(src_img);title('原图像');  
[a,b]=size(src_img);
dst_img=zeros(a,b);

var=0.05;
%noise_type='salt & pepper';
noise_type='speckle';
%noise_type='gaussian';
img_noise=imnoise(src_img,noise_type,var);

%%%%%%%%进行叠加平均%%%%%%%%%%%%%%%%%%%%
subplot(232)
imshow(img_noise);title('噪声图像');
for n=1:200
    img_noise=imnoise(src_img,noise_type,var);
    dst_img=dst_img+double(img_noise);%求和
end;
dst_img=dst_img/n;%取平均
%求取叠加平均除燥后的psnr结果
pnsr_averge = myPsnr(dst_img,src_img)
%求取加燥后的psnr结果
pnsr_noise = myPsnr(img_noise,src_img)

dst_img=uint8(dst_img);
subplot(233)
imshow(dst_img);title('叠加平均');

%%%%%%%%进行小波去噪%%%%%%%%%%%%%%%%%%%%
% 获取输入参数    
w    = 'db6';%小波类型    
n    = 3;%分解层数     
  
sorh1 = 'hard';%硬阈值    
sorh2 = 'soft';%软阈值    
sorh3 = 'half';%半阈值     
% 对图像进行小波分解    
[c,l] = wavedec2(img_noise,n,w);    
    
  
%求取阈值    
N = numel(src_img);    
[chd1,cvd1,cdd1] = detcoef2('all',c,l,1);     
cdd1=cdd1(:)';    
sigma = median(abs(cdd1))/0.6745;%提取细节系数求中值并除以0.6745    
thr = sigma*sqrt(2*log(N))/1.414; %对阈值做了改进  
  
  
% 对小波系数全局阈值处理    
cxchard = c;% 保留近似系数    
cxcsoft = c;% 保留近似系数    
cxchalf = c;% 保留近似系数    
justdet = prod(l(1,:))+1:length(c);%截取细节系数（不处理近似系数）    
    
% 阈值处理细节系数    
cxchard(justdet) = myWthresh(cxchard(justdet),sorh1,thr);%硬阈值去噪    
cxcsoft(justdet) = myWthresh(cxcsoft(justdet),sorh2,thr);%软阈值去噪    
cxchalf(justdet) = myWthresh(cxchalf(justdet),sorh3,thr);%软阈值去噪    
  
%小波重建    
xchard = waverec2(cxchard,l,w);    
xcsoft = waverec2(cxcsoft,l,w);    
xchalf = waverec2(cxchalf,l,w);     
  
subplot(234)   
imshow(uint8(xchard));title('硬阈值去噪图像');      
pnsr_hard = myPsnr(xchard,src_img)

subplot(235)   
imshow(uint8(xcsoft));title('软阈值去噪图像');     
pnsr_soft = myPsnr(xcsoft,src_img)

subplot(236)   
imshow(uint8(xchalf));title('半阈值去噪图像');   
pnsr_half = myPsnr(xchalf,src_img)

pnsr程序

function pnsr_result = myPsnr(img_ref,img_den)
[a,b]=size(img_ref);
XX2=double(img_ref) - double(img_den);
mse_value2 = sum(sum( XX2.^2 ))/(a*b);
pnsr_result = 10*log10( 255*255 / mse_value2 );