DIP关键算法-去噪算法评价标准

　　去噪算法是图像处理里面一个基础的问题，对于去噪算法的好坏一个直接的判断就是视觉上的感知，但是这只能作为一个定性的判断，在实际的应用中，我们还需要进行定量分析，进而做出最终的判断。在实际中往往采用原始图像Iorigin和去噪图像Idenoise之间的残差图像Iresidual=Iorigin−Idenoise作为定量分析的手段，当残差图像没有明显的结构的时候，即残差图像杂乱无章，类似于高斯白噪声，我们就可以认为去噪算法效果比较好，反之，则认为效果比较差。如下是一个非局部均值去噪算法的结果分析。

　　


　　

代码如下

clear
clc
clf
colormap(gray)

% create example image
ima=100*ones(100);
ima(50:100,:)=50;
ima(:,50:100)=2*ima(:,50:100);
fs=fspecial('average');
ima=imfilter(ima,fs,'symmetric');

% add some noise
sigma=10;
rima=ima+sigma*randn(size(ima));

% show it
imagesc(rima)
drawnow

% denoise it
fima=NLmeansfilter(ima,5,2,sigma);

% show results
clf
subplot(2,2,1),imagesc(ima),title('original');
subplot(2,2,2),imagesc(rima),title('noisy');
subplot(2,2,3),imagesc(fima),title('filtered');
subplot(2,2,4),imagesc(rima-fima),title('residuals');

其中的算法源码如下

function [output]=NLmeansfilter(input,t,f,h)

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
%  input: image to be filtered
%  t: radio of search window
%  f: radio of similarity window
%  h: degree of filtering
%
%  Implementation of the Non local filter proposed for A. Buades, B. Coll and J.M. Morel in
%  "A non-local algorithm for image denoising"
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% Size of the image
[m n]=size(input);

% Memory for the output
Output=zeros(m,n);

% Replicate the boundaries of the input image
input2 = padarray(input,[f f],'symmetric');

% Used kernel
kernel = make_kernel(f);
kernel = kernel / sum(sum(kernel));

h=h*h;

for i=1:m
for j=1:n

i1 = i+ f;
j1 = j+ f;

W1= input2(i1-f:i1+f , j1-f:j1+f);

wmax=0;
average=0;
sweight=0;

rmin = max(i1-t,f+1);
rmax = min(i1+t,m+f);
smin = max(j1-t,f+1);
smax = min(j1+t,n+f);

for r=rmin:1:rmax
for s=smin:1:smax

if(r==i1 && s==j1) continue; end;

W2= input2(r-f:r+f , s-f:s+f);

d = sum(sum(kernel.*(W1-W2).*(W1-W2)));

w=exp(-d/h);

if w>wmax
wmax=w;
end

sweight = sweight + w;
average = average + w*input2(r,s);
end
end

average = average + wmax*input2(i1,j1);
sweight = sweight + wmax;

if sweight > 0
output(i,j) = average / sweight;
else
output(i,j) = input(i,j);
end
end
end

function [kernel] = make_kernel(f)

kernel=zeros(2*f+1,2*f+1);
for d=1:f
value= 1 / (2*d+1)^2 ;
for i=-d:d
for j=-d:d
kernel(f+1-i,f+1-j)= kernel(f+1-i,f+1-j) + value ;
end
end
end
kernel = kernel ./ f;