图像复原MATLAB实现

前言：本篇博客先介绍滤波器滤除噪声，再介绍滤波器复原，侧重于程序的实现。

一：三种常见的噪声

 

 

 

 二：空间域滤波

空间域滤波复原是在已知噪声模型的基础上，对噪声的空间域进行滤波。
空间域滤波复原方法主要包括：

　　均值滤波器

　　　　算术均值滤波器

　　　　几何均值滤波器
　　　　谐波均值滤波器
　　　　逆谐波均值滤波器
　　顺序统计滤波器
　　　　中值滤波器
　　　　最大值/最小值滤波器

2.1算数均值滤波器

img=imread('D:/picture/ZiXia.jpg');
img=rgb2gray(img);
figure,imshow(img);//原图
img_noise=double(imnoise(img,'gaussian',0.06));
figure,imshow(img_noise,[]);//含有高斯噪声的图
img_mean=imfilter(img_noise,fspecial('average',3));//滤波后的图
figure;imshow(img_mean,[]);

2.2几何均值滤波器

img=imread('cameraman.tif');
img=rgb2gray(img);
figure,imshow(img);
img_noise=double(imnoise(img,'gaussian',0.06));
figure,imshow(img_noise,[]);
img_mean=exp(imfilter(log(img_noise+1),fspecial('average',3)));
figure;imshow(img_mean,[]);

 

 

 2.3谐波均值滤波器

2.4逆谐波均值滤波器

 

采用逆谐波均值滤波器对附加胡椒噪声图像进行滤波的matlab程序如下：

img=imread('cameraman.tif');  figure,imshow(img);
[M,N]=size(img);R=imnoise2('salt & pepper',M,N,0.1,0);
img_noise=img;img_noise(R==0)=0;
img_noise=double(img_noise); figure,imshow(img_noise,[]);
Q=1.5;
img_mean=imfilter(img_noise.^(Q+1),fspecial('average',3))./imfilter(img_noise.^Q,fspecial('average',3));
figure;imshow(img_mean,[]);

采用逆谐波均值滤波器对附加盐噪声图像进行滤波的matlab程序如下：

img=imread('csboard.tif');figure,imshow(img);
[M,N]=size(img);R=imnoise2('salt & pepper',M,N,0,0.1);
img_noise=img;img_noise(R==1)=255;
img_noise=double(img_noise); figure,imshow(img_noise,[]);
Q=-1.5;
img_mean=imfilter(img_noise.^(Q+1),fspecial('average',3))./imfilter(img_noise.^Q,fspecial('average',3));
figure;imshow(img_mean,[]);

 2.5中值滤波器

 

 

img=imread('cameraman.tif');
img_noise=double(imnoise(img,'salt & pepper',0.06));
img_mean=imfilter(img_noise,fspecial('average',5));
img_median=medfilt2(img_noise);%一次中值滤波
img_median2=medfilt2(img_median);%二次中值滤波

2.6最大值，最小值滤波器

 

 

利用最大值滤波器消除胡椒噪声污染图像的matlab程序如下。

img=imread('csboard.tif');
[M,N]=size(img);
R=imnoise2('salt & pepper',M,N,0.1,0);
img_noise=img;
img_noise(R==0)=0;
img_noise=double(img_noise);
imwrite(uint8(img_noise),'csbord_pepper.jpg');
img_max=imdilate(img_noise,ones(3,3));
imwrite(uint8(img_max),'cameraman_saltpepper_max.jpg');

利用最小值滤波器消除盐噪声污染图像的matlab程序如下。

img=imread('csboard.tif');
[M,N]=size(img);
R=imnoise2('salt & pepper',M,N,0,0.1);
img_noise=img;
img_noise(R==1)=255;
img_noise=double(img_noise);
imwrite(uint8(img_noise),'csbord_salt.jpg');
img_min=imerode(img_noise,ones(3,3));
imwrite(uint8(img_min),'cameraman_saltpepper_min.jpg');

 2.7带阻滤波器

 

 

I=imread('pout.tif');
[m,n]=size(I);
J=I;
for i=1:m
    for j=1:n
        J(i,j)=I(i,j)+20*sin(20*i)+20*sin(20*j);%增加周期性噪声
    end
end
IF=fftshift(fft2(I));
JF=fftshift(fft2(J));
IF=log(1+abs(IF));
JF=log(1+abs(JF));

subplot(221)%显示频谱
imshow(IF,[])
subplot(222)
imshow(JF,[])
%高斯带阻滤波器构造
fbrf=ones(m,n);
for i=1:m
    for j=1:n
        fbrf(i,j)=1-exp(-0.5*(((i-m/2)^2+(j-n/2)^2-55^2)/(sqrt(i.^2+j.^2)*5))^2);%20为带阻中心，5为带宽
    end
end
 H=fbrf;
%频率域滤波
f=fftshift(fft2(J));
out=f.*H;%频率域滤波
out=ifft2(ifftshift(out));
out=abs(out);
out=out/max(out(:));%归一化【0,1】
subplot(223)
imshow(out,[]);
subplot(224)
imshow(J,[])

2.8带通滤波器

2.9陷波滤波器

 

clear;
close all;
src = im2double(imread('D:/picture/ZiXia.jpg'));
src = rgb2gray(src);
subplot(221)
imshow(src);
title('原始图像');
[w h] = size(src);
srcf = fft2(src);
srcf = fftshift(srcf);
subplot(222)
imshow(srcf);
%  低通滤波
% flt = zeros(size(src));
% rx1 = w/2;
% ry1 = h/2;
% r = min(w,h)/3;
% for i = 1:w
%     for j = 1:h
%         if(rx1-i)^2 +(ry1 - j)^2 <= r*r
%             flt(i,j) = 1;
%         end
%     end
% end
% 陷波滤波
flt = ones(size(src));
r = min(w,h)/12;
rx1 = r
ry1 =h/2
for i = 1:w
    for j = 1:h
        if(rx1-i)^2 +(ry1 - j)^2 <= r*r
            flt(i,j) = 0;
        end
        if(w-rx1-i)^2 +(h-ry1 - j)^2 <= r*r
            flt(i,j) = 0;
        end
    end
end
subplot(223)
imshow(flt);
title('滤波器图像');
dfimg = srcf.*flt;
dfimg = ifftshift(dfimg);
dimg = ifft2(dfimg,'symmetric');
subplot(224)
imshow(dimg):title('滤波后');

 

 

 2.9逆滤波

image_o=imread('D:/picture/lenagray.jpg');
subplot(1,3,1);
imshow(image_o);
title('原图像');
%频率域退化图像，退化函数H(u,v)=exp(-0.0025*( (u-M/2).^2+(v-N/2).^2).^(5/6) )
%傅里叶变换
f=im2double(image_o);
F=fft2(f);
F=fftshift(F);
%执行退化
[M,N]=size(F);
[u,v]=meshgrid(1:M,1:N);%生成二维坐标系
H=exp(-0.0025* ( (u-M/2).^2+(v-N/2).^2).^(5/6) );
F=F.*H;
%傅里叶反变换
X=ifftshift(F);
x=ifft2(X);
x=uint8(abs(x)*256);
subplot(1,3,2);
imshow(x);
%
title('退化图像');

image_d=x;
%直接逆滤波图像复原

ff=im2double(image_d);%将图像灰度值归一化到0-1之间

% 傅里叶变换
f_Id=fft2(ff);
f_Id=fftshift(f_Id);
fH_Id=f_Id;
[M,N]=size(fH_Id);
% 逆滤波
threshold=78;
if threshold>M/2
        %全滤波
        fH_Id=fH_Id./(H+eps);
else
        %对一定半径范围内进行滤波
        for i=1:M
            for j=1:N
                if sqrt((i-M/2).^2+(j-N/2).^2)<threshold
                    fH_Id(i,j)=fH_Id(i,j)./(H(i,j)+eps);
                end
            end
        end
end

% 执行傅立叶逆变换
fH_Id1=ifftshift(fH_Id);
f_new=ifft2(fH_Id1);
f_new=uint8(abs(f_new)*255);
subplot(1,3,3);
imshow(f_new);
title('滤波半径=78的逆滤波复原图像');

2.10维纳滤波

 

 

 直接截图了，没敲