数字图像处理实验四

5.4 实验四 图像的频率域滤波(综合性)

一、实验目的

1.掌握数字图像傅里叶变换;
2.掌握基于频率域的图像滤波算法;
3.能够基于MATLAB实现基于频率域的图像滤波;
4.能够根据实验结果分析各种算法的特点以及应用场合,培养处理实际图像的能力。

二、实验要求

1.实验课前需要写预习实验报告,内容为本次实验要求中的所有程序清单。
2.实验课对预习报告中的编程代码进行上机调试,完成实验指导书中全部实验要求内容。
3.实验课后写出实验报告。报告要求有实验目的,实验内容与步骤,调试完成的准确编程代码,实验小结,回答问题。

三、实验内容(每一个内容编写一个*.m文件)

1.用MATLAB函数对图像(house.tif)完成傅立叶变换;
(1)读入,在图形窗口显示输入图像;
代码:

clc;
clf;
I = imread("C:\Users\25023\Desktop\图像识别文件\实验四用到的图像\house.tif");
figure(1);
subplot(1,1,1)
imshow(I);
title("原图")

结果:
在这里插入图片描述

(2)完成傅立叶变换,显示变换后的频谱图;
代码:

clc;
clf;
I = imread("C:\Users\25023\Desktop\图像识别文件\实验四用到的图像\house.tif");
figure(1);
 
subplot(1,2,1)
imshow(I);
I = im2double(I);
title("原图")
%f=rgb2gray(img);    %对于RGB图像必须做的一步,也可以用im2double函数
 
 
subplot(1,2,2);
F=fft2(I);          %傅里叶变换
imshow(uint8(abs(F)));   %abs取复数矩阵的模
title('傅里叶变换后的频谱图');

结果:
在这里插入图片描述

(3)对频谱图采用对数变换,将变换后的图像在图形窗口显示;
代码:

clc;
clf;
I = imread("C:\Users\25023\Desktop\图像识别文件\实验四用到的图像\house.tif");
figure(1);
 
subplot(2,2,1)
imshow(I);
I = im2double(I);
title("原图")
%f=rgb2gray(img);    %对于RGB图像必须做的一步,也可以用im2double函数
 
 
subplot(2,2,2);
F=fft2(I);          %傅里叶变换
imshow(uint8(abs(F)));   %abs取复数矩阵的模
title('傅里叶变换后的频谱图');
 
subplot(2,2,3);
F1=log(abs(F)+1);   %取模并进行缩放
imshow(F1,[]);
title('傅里叶对数变换频谱图');  

结果:
在这里插入图片描述

(4)对变换后的频谱进行频移,将零频率部分移到图形的中心位置,并在图形窗口显示;
代码:

clc;
clf;
I = imread("C:\Users\25023\Desktop\图像识别文件\实验四用到的图像\house.tif");
figure(1);
 
subplot(2,3,1)
imshow(I);
I = im2double(I);
title("原图")
%f=rgb2gray(img);    %对于RGB图像必须做的一步,也可以用im2double函数
 
 
subplot(2,3,2);
F=fft2(I);          %傅里叶变换
imshow(abs(F));   %abs取复数矩阵的模
title('傅里叶变换后的频谱图');
 
subplot(2,3,3);
F1=log(abs(F)+1);   %取模并进行缩放
imshow(F1,[]);
title('傅里叶对数变换频谱图');  
 
 
subplot(2,3,4);
Fs=fftshift(F1);      %将频谱图中零频率成分移动至频谱图中心
imshow(abs(Fs),[]);
title('频移后的频谱图');

结果:
在这里插入图片描述

(5)对频谱图采用对数变换,将变换后的图像在图形窗口显示;
代码:

clc;
clf;
I = imread("C:\Users\25023\Desktop\图像识别文件\实验四用到的图像\house.tif");
figure(1);
 
subplot(2,3,1)
imshow(I);
I = im2double(I);
title("原图")
%f=rgb2gray(img);    %对于RGB图像必须做的一步,也可以用im2double函数
 
 
subplot(2,3,2);
F=fft2(I);          %傅里叶变换
imshow(abs(F));   %abs取复数矩阵的模
title('傅里叶变换后的频谱图');
 
subplot(2,3,3);
F1=log(abs(F)+1);   %取模并进行缩放
imshow(F1,[]);
title('傅里叶对数变换频谱图');  
 
 
subplot(2,3,4);
Fs=fftshift(F1);      %将频谱图中零频率成分移动至频谱图中心
imshow(abs(Fs),[]);
title('频移后的频谱图');
 
subplot(2,3,5);
Fs=fftshift(F);      %将频谱图中零频率成分移动至频谱图中心
S=log(abs(Fs)+1);    %取模并进行缩放
imshow(S,[]);
title('频移后的对数变换频谱图');

结果:
在这里插入图片描述

(6)获取傅里叶变换后的相位谱图,并在图形窗口显示;
代码:

clc;
clf;
I = imread("C:\Users\25023\Desktop\图像识别文件\实验四用到的图像\house.tif");
figure(1);
 
subplot(2,3,1)
imshow(I);
I = im2double(I);
title("原图")
%f=rgb2gray(img);    %对于RGB图像必须做的一步,也可以用im2double函数
 
 
subplot(2,3,2);
F=fft2(I);          %傅里叶变换
phase=angle(F)*180/pi;     %图像相位谱
imshow(phase,[]);
title('傅里叶变换后的相位图');

结果:
在这里插入图片描述

(7)对相位谱图采用对数变换,将变换后的图像在图形窗口显示;
代码:

clc;
clf;
I = imread("C:\Users\25023\Desktop\图像识别文件\实验四用到的图像\house.tif");
figure(1);
 
subplot(1,3,1)
imshow(I);
I = im2double(I);
title("原图")
%f=rgb2gray(img);    %对于RGB图像必须做的一步,也可以用im2double函数
 
 
subplot(1,3,2);
F=fft2(I);          %傅里叶变换
phase=angle(F)*180/pi;     %图像相位谱
imshow(phase,[]);
title('傅里叶变换后的相位图');
 
subplot(1,3,3);
phase_2=log(angle(F)*180/pi);     %图像相位谱
imshow(phase_2,[]);
title('傅里叶变换后的相位图的对数相位图');

结果:
在这里插入图片描述
(8)通过傅里叶反变换得到输入的图像,显示反变换得到的图像。
代码:

clc;
clf;
I = imread("C:\Users\25023\Desktop\图像识别文件\实验四用到的图像\house.tif");
figure(1);
 
subplot(2,3,1)
imshow(I);
I = im2double(I);
title("原图")
%f=rgb2gray(img);    %对于RGB图像必须做的一步,也可以用im2double函数
 
 
subplot(2,3,2);
F=fft2(I);          %傅里叶变换
imshow(abs(F));   %abs取复数矩阵的模
title('傅里叶变换后的频谱图');
 
subplot(2,3,3);
F1=log(abs(F)+1);   %取模并进行缩放
imshow(F1,[]);
title('傅里叶对数变换频谱图');  
 
 
subplot(2,3,4);
Fs=fftshift(F);      %将频谱图中零频率成分移动至频谱图中心
imshow(abs(Fs),[]);
title('频移后的频谱图');
 
 
subplot(2,3,5);
fr=real(ifft2(ifftshift(Fs)));  %频率域反变换到空间域,并取实部
ret=im2uint8(mat2gray(fr));    %更改图像类型
imshow(ret);
title('逆傅里叶变换');

结果:
在这里插入图片描述

2.用MATLAB生成如下的两幅图像(图a:黑白线竖线8条,图b: 黑白线横线32条),图像的大小为,对这两幅图像进行傅里叶变换,并分别进行频移和对数变换,在一个图形窗口显示:原图像、傅里叶变换后的频谱图和相位图、频移后的频谱图及相位图。
在这里插入图片描述

图a:黑白线竖线8条 图b: 黑白线横线32条
代码:

clc,
clear
lines = 8; %黑白线对数
m = 16;       %黑或白线宽度(像素),黑白宽度一致
n = 256;     %图像高度
S = zeros(1,2*m);
S(1,1:m) = ones(1,m);
A = repmat(S,n,lines);
subplot(3,4,1),imshow(A),title('图一');
 
 
figure(1)
F=fft2(A);
subplot(2,4,2),imshow(abs(F)),title('傅里叶变换后的频谱');
T=log(F+1);
subplot(2,4,3),imshow(abs(T)),title('对数变换后');
S = fftshift(F);
subplot(2,4,4),imshow(abs(S)),title('频谱中心');
subplot(2,4,5),imshow(abs(log(S+1))),title('对数变换频谱中心');
subplot(2,4,6),imshow((angle(F)*180/pi),[]),title('相位谱图1');
subplot(2,4,7),imshow((log(angle(F)*180/pi)+1),[]),title('相位谱图2');
 
figure(2)
lines = 32; %黑白线对数
m = 4;       %黑或白线宽度(像素),黑白宽度一致
n = 256;     %图像高度
S = zeros(2*m,1);
S(1:m,1) = ones(1,m);
B = repmat(S,lines,n);
subplot(2,4,1),imshow(B),title('图二');
 
F=fft2(B);
subplot(2,4,2),imshow(abs(F)),title('傅里叶变换后的频谱');
T=log(F+1);
subplot(2,4,3),imshow(abs(T)),title('对数变换后');
S = fftshift(F);
subplot(2,4,4),imshow(abs(S)),title('频谱中心');
subplot(2,4,5),imshow(abs(log(S+1))),title('对数变换频谱中心');
subplot(2,4,6),imshow((angle(F)*180/pi),[]),title('相位谱图1');
subplot(2,4,7),imshow((log(angle(F)*180/pi)+1),[]),title('相位谱图2');

结果:
在这里插入图片描述

3.比较傅里叶变换后图像的频谱图和相位图所起的作用。
(1)读入并显示两幅不同的彩色人体图像(man.png和woman.png),大小相同;

(2)分别对上述两幅图像进行傅里叶变换并频移,得到他们的频谱图和相位图;

(3)根据第一个人员图像的频谱图和第二个人员图像的相位图恢复图像;

(4)根据第一个人员图像的相位图和第二个人员图像的频谱图恢复图像;

(5)比较分析图像与原图的差别。

clc;
clf;
I=imread('C:\Users\25023\Desktop\图像识别文件\实验四用到的图像\man.png');
J=imread('C:\Users\25023\Desktop\图像识别文件\实验四用到的图像\woman.png');
I = rgb2gray(I);
J = rgb2gray(J);%转灰度图像
[M,N]=size(I);
J=imresize(J,[M,N]);%限制图像大小,相位交换需要两个图像大小一致
 
IF1 = fft2(I);%图像进行傅里叶变换
JF1 = fft2(J);%图像进行傅里叶变换
 
IF2 = fftshift(abs(IF1));%图像频谱中心化
JF2 = fftshift(abs(JF1));%图像频谱中心化
 
IABS = abs(IF1);
IAG = angle(IF1);%幅度谱,相位谱
JABS = abs(JF1);
JAG = angle(JF1);
 
I2=fftshift(abs(IABS)); %log(abs(I2)+1)等于log(abs(IF2)+1)
J2=fftshift(abs(JABS));
figure(1)
 
subplot(2,4,1);
imshow(I);
title('原始图像1');
 
subplot(2,4,2);
imshow(log(IABS+1),[]);
title('图像1的幅度谱');
 
subplot(2,4,3);
imshow(IAG,[]);
title('图像1的相位谱');
 
subplot(2,4,5);
imshow(J);title('原始图像2');
 
subplot(2,4,6);
imshow(log(JABS),[]);
title('图像2的幅度谱');
 
subplot(2,4,7);
imshow(JAG,[]);
title('图像2的相位谱');
 
subplot(2,4,4);
imshow(log(abs(I2)+1),[]);
title('图像1的频谱中心化');
 
subplot(2,4,8);
imshow(log(abs(J2)+1),[]);
title('图像2的频谱中心化');
 
II=ifft2(IF1);
JJ=ifft2(JF1);
figure(2)
 
subplot(1,2,1);
imshow(abs(II),[]);
title('图像1逆变换');
 
subplot(1,2,2);
imshow(abs(JJ),[]);
title('图像2逆变换');
 
%相位交换,注意i,两图像一定要大小一样
ii=IABS.*cos(JAG)+IABS.*sin(JAG).*1i;
jj=JABS.*cos(IAG)+JABS.*sin(IAG).*1i;
iI=abs(ifft2(ii));
jJ=abs(ifft2(jj));
figure(3)
 
subplot(1,2,1);
imshow(iI,[]);
title('图像1 幅度谱 与 图像2相位谱');
 
subplot(1,2,2);
imshow(jJ,[]);
title('图像2 幅度谱 与 图像1相位谱')

结果:
在这里插入图片描述

4.低通滤波和高通滤波的效果比较
(1)读入并显示图像Fig0313(a).tif;
(2)分别对图像进行理想低通滤波和高斯低通滤波,D0=20,60,将以上结果在一个图形窗口显示;
(3)分别对图像进行理想高通滤波和高斯高通滤波,D0=20,60,将以上结果在一个图形窗口显示;
(4)比较分析滤波结果的区别。
代码:

clc;
clf;
I  = imread("C:\Users\25023\Desktop\图像识别文件\实验四用到的图像\Fig0313(a).tif");
figure(1);
subplot(1,5,1),imshow(I);
title('原图像');
s = lixiang_low(I);
subplot(1,5,2),imshow(s);
title('理想低通滤波所得图像 D0=20'); 
 
d = lixiang_high(I);
subplot(1,5,3),imshow(d);
title('理想高通滤波所得图像 D0=20'); 
 
Ia = gaosi_low(I);
subplot(1,5,4),imshow(uint8(Ia));
title('高斯低通滤波处理 D0=60');
 
img_lpf = gaosi_high(I);
subplot(1,5,5);
imshow(img_lpf);title('高斯高通滤波D0=60');
 
function result = lixiang_low(I)
s=fftshift(fft2(I));
[a,b]=size(s);
a0=round(a/2);
b0=round(b/2);
d=20;
for i=1:a 
    for j=1:b 
        distance=sqrt((i-a0)^2+(j-b0)^2);
        if distance<=d
            h=1;
        else
            h=0;
        end
        s(i,j)=h*s(i,j);
    end
end
s=uint8(real(ifft2(ifftshift(s))));
result = s;
end
 
function result = lixiang_high(I)
s=fftshift(fft2(I));
[a,b]=size(s);
a0=round(a/2);
b0=round(b/2);
d=20;
for i=1:a 
    for j=1:b 
        distance=sqrt((i-a0)^2+(j-b0)^2);
        if distance>=d
            h=1;
        else
            h=0;
        end
        s(i,j)=h*s(i,j);
    end
end
s=uint8(real(ifft2(ifftshift(s))));
result = s;
end
 
function result = gaosi_low(IA)
[f1,f2]=freqspace(size(IA),'meshgrid');
D=60/size(IA,1);
r=f1.^2+f2.^2;
Hd=ones(size(IA));
for i=1:size(IA,1)
    for j=1:size(IA,2)
        t=r(i,j)/(D*D);
        Hd(i,j)=exp(-t);
    end
end
Y=fft2(double(IA));
Y=fftshift(Y);
Ya=Y.*Hd;
Ya=ifftshift(Ya);
Ia=real(ifft2(Ya));
result = Ia;
end
 
 
function result = gaosi_high(image)
d0=60;  %阈值
[M ,N]=size(image);
 
img_f = fft2(double(image));%傅里叶变换得到频谱
img_f=fftshift(img_f);  %移到中间
 
m_mid=floor(M/2);%中心点坐标
n_mid=floor(N/2);  
 
h = zeros(M,N);%高斯低通滤波器构造
for i = 1:M
    for j = 1:N
        d = ((i-m_mid)^2+(j-n_mid)^2);
        h(i,j) = 1-exp(-(d)/(2*(d0^2)));      
    end
end
img_lpf = h.*img_f;
img_lpf=ifftshift(img_lpf);    %中心平移回原来状态
img_lpf=uint8(real(ifft2(img_lpf)));  %反傅里叶变换,取实数部分
result = img_lpf;
end

结果:

在这里插入图片描述

5.(选做)空间域滤波和频率域滤波的比较
(1)读入图像house.tif,对其进行傅里叶变换,显示原图像和傅里叶变换后的频谱图像;
(2)用MATLAB函数fspecial()生成空间滤波器(sobel,laplacian),对图像进行滤波,显示滤波后的图像;
(3)用MATLAB函数将空间滤波器转化为频率域滤波器,对图像进行频率域的滤波;
(4)比较空间滤波和频率域滤波速度和滤波后图像的区别。

  • 4
    点赞
  • 23
    收藏
    觉得还不错? 一键收藏
  • 1
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论 1
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值