实验四 图像频域平滑与锐化

实验内容

1、图像频域平滑（去噪）：使用自生成图像（包含白色区域，黑色区域，并且部分区域添加椒盐噪声），然后进行傅里叶变换，并且分别使用理想低通滤波器、巴特沃斯低通滤波器、指数低通滤波器和梯形低通滤波器（至少使用两种低通滤波器），显示滤波前后的频域能量分布图，空间图像。分析不同滤波器对噪声、边缘的处理效果及其优缺点。

data = rand(50,50,3);%自生成一个图像矩阵

data=zeros(50);

for i=1:30

    for j =1:30

        data(i,j)=0;    %randi([180, 255]);

    end

end

for i=40:50

    for j=40:50

        data(i,j)=1;%白

    end

end

K =imnoise(data,'salt & pepper',0.04); %加椒盐噪声

subplot(2,3,1);

imshow(K);

title('自生成图像');

 

s=fftshift(fft2(K));%先傅里叶变换，然后将傅里叶变化零频率搬移到频谱中间

subplot(2,3,2);

imshow(abs(s),[]);%abs取复数的幅值，[]自动调整数据的范围以便于显示

title('图像傅里叶变换所得频谱');

[a,b]=size(s);

a0=round(a/2);%确定傅里叶变化原点，并且数据取整

b0=round(b/2);

d0=15;%截止频率

d1=30;

n=2;

s1=zeros(a,b);%理想低通滤波

s2=zeros(a,b);%巴特沃斯低通滤波

s3=zeros(a,b);%指数低通滤波

for i=1:a

    for j=1:b

        d=sqrt((i-a0)^2+(j-b0)^2);

        if d<=d0

            h1=1;

        else

            h1=0;

        end

        h2=1/(1+0.414*((d/d0)^(2*n)));

        h3=exp(-0.347*(d/d0)^n);

        if d<=d0

            h4=1;

        elseif d<=d1

            h4=(d-d1)/(d0-d1);

        else

            h4=0;

        end

        s1(i,j)=h1*s(i,j);

        s2(i,j)=h2*s(i,j);

        s3(i,j)=h3*s(i,j);

        s4(i,j)=h4*s(i,j);

    end

end

s1=real(ifft2(ifftshift(s1)));%傅里叶逆移频，由于之前做过fftshift，然后傅里叶反变换，然后提取实部

subplot(2,3,3);

imshow(s1);

title('理想低通滤波所得图像');

s2=real(ifft2(ifftshift(s2)));

subplot(2,3,4);

imshow(s2);

title('巴特沃斯低通滤波所得图像');

s3=real(ifft2(ifftshift(s3)));

subplot(2,3,5);

imshow(s3);

title('指数低通滤波所得图像');

2、图像频域平滑（锐化）：选择一幅图像，例如rice.png，分别使用理想高通滤波器、巴特沃斯高通滤波器、指数高通滤波器和梯形高通滤波器（至少使用两种高通滤波器），显示滤波前后的频域能量分布图，空间图像。分析不同滤波器处理效果及其优缺点。

I=imread('L:\图像处理\实验文件\rice.png');

subplot(3,4,1);

imshow(I);

title('原图像');

s=fftshift(fft2(I));%先傅里叶变换，然后将傅里叶变化零频率搬移到频谱中间

subplot(3,4,2);

imshow(log(abs(s)),[]);

title('原图频谱');

[a,b]=size(s);

a0=round(a/2);%确定傅里叶变化原点，并且数据取整

b0=round(b/2);

d0=15;%截止频率

d1=30;

n=2;

s1=zeros(a,b);%理想高通滤波

s2=zeros(a,b);%巴特沃斯高通滤波

s3=zeros(a,b);%指数高通滤波

for i=1:a

    for j=1:b

        d=sqrt((i-a0)^2+(j-b0)^2);

        if d<=d0

            h1=0;

        else

            h1=1;

        end

        h2=1/(1+(d0/d)^(2*n));

        h3=exp(-(d0/d)^n);

        if d<=d0

            h4=0;

        elseif d<=d1

            h4=(d-d1)/(d0-d1);

        else

            h4=1;

        end

        s1(i,j)=h1*s(i,j);

        s2(i,j)=h2*s(i,j);

        s3(i,j)=h3*s(i,j);

        s4(i,j)=h4*s(i,j);

    end

end

s1=uint8(real(ifft2(ifftshift(s1))));%傅里叶逆移频，由于之前做过fftshift，然后傅里叶反变换，然后提取实部

subplot(3,4,3);

imshow(s1);

title('理想高通滤波所得图像');

s1new=fftshift(fft2(s1));

subplot(3,4,4);

imshow(log(abs(s1new)),[]);

title('理想高通滤波新频谱');

s2=uint8(real(ifft2(ifftshift(s2))));

subplot(3,4,5);

imshow(s2);

title('巴特沃斯高通滤波所得图像');

s2new=fftshift(fft2(s2));

subplot(3,4,6);

imshow(log(abs(s2new)),[]);

title('巴特沃斯高通滤波新频谱');

s3=uint8(real(ifft2(ifftshift(s3))));

subplot(3,4,7);

imshow(s3);

title('指数高通滤波所得图像');

s3new=fftshift(fft2(s3));

subplot(3,4,8);

imshow(log(abs(s3new)),[]);

title('指数高通滤波新频谱');