【数字信号处理】实验三 离散傅里叶变换

一、 实验目的

应用MATLAB数字信号处理工具箱处理离散信号（DFT），对各种典型序列做傅立叶变换，验证N点DFT的物理意义，严整频域采样与时域采样的对偶性及快速卷积算法

二、 实验内容

（1）基本序列的离散傅立叶变换计算

clear;close all
N=16;N1=8;n=0:N-1;
k=0:N1-1;x1n=exp(j*pi*n/8);
X1k=fft(x1n,N);Xk1=fft(x1n,N1);
x2n=cos(pi*n/8);X2k=fft(x2n,N);
Xk2=fft(x2n,N1);x3n=sin(pi*n/8);
X3k=fft(x3n,N);Xk3=fft(x3n,N1);
subplot(2,3,1),stem(n,abs(X1k),'.')
subplot(2,3,2),stem(n,abs(X2k),'.')
subplot(2,3,3),stem(n,abs(X3k),'.')
subplot(2,3,4),stem(k,abs(Xk1),'.')
subplot(2,3,5),stem(k,abs(Xk2),'.')
subplot(2,3,6),stem(k,abs(Xk3),'.')

（2）验证N点DFT的物理意义

clear;close all  %关闭所有图像并清除已有数据
N1=8;N2=16;n=0:N1-1;w=2*pi* (0:2047)/2048;
Xw=(1-exp(-j*4*w))./(1-exp(-j*w));
subplot(3,2,1);plot(w/pi,abs(Xw))title('x(n)的相频曲线');
xlabel('ω/π');ylabel('幅度');
subplot(3,2,2);plot(w/pi,angle(Xw))title('x(n)的相频曲线');
axis([0,2,-pi,pi]);line([0,2],[0,0])xlabel('ω/π');ylabel('相位(rad)');
xn=[(n>=0)&(n<4)];
X1k=fft(xn,N1);
X2k=fft(xn,N2);
k1=0:N1-1;
subplot(3,2,3);stem(k1,abs(X1k),'.')title('N1 点 DFT[x(n)]=X1(k)')xlabel('k (ω=2πk/N1)');ylabel('|X1(k)|');
hold on
plot(N1/2*w/pi,abs(Xw))
subplot(3,2,4);stem(k1,angle(X1k),'.')title('X1(k)的相位');
axis([0,N1,-pi,pi]);
line([0,N1],[0,0])
xlabel('k(ω=2πk/N1)');ylabel('相位(rad)');
hold on   %保持原图，并支持叠加绘图
plot(N1/2*w/pi,angle(Xw))
k2=0:N2-1;
subplot(3,2,5);
stem(k2,abs(X2k),'.')title('N2 点 DFT[x(n)]=X2(k)')xlabel('k(ω=2πk/N2)'):ylabel('|X2(k)|');
hold on %保持原图，并支持叠加绘图
plot(N2/2*w/pi,abs(Xw))
subplot(3,2,6);
stem(k2,angle(X2k),'.')title('N2 点 DFT[x(n)]=X2(k)')
hold on  %保持原图，并支持叠加绘图
plot(N2/2*w/pi,angle(Xw))

（3）快速卷积

clear;close all
xn=input('请输入 x(n)序列：xn= ');
hn=input('请输入 h(n)长度:hn= ');
M=length(xn);
N=length(hn);
nx=1:M; nh=1:N;
L=pow2(nextpow2(M+N-1)); %取L为≥且接近(M+N-1)的2正次幂
tic,  %快速卷积的开始
Xk=fft(xn,L);Hk=fft(hn,L);
Yk=Xk.*Hk;yn=ifft(Yk,L);toc  %快速卷积的结束
subplot(2,2,1),stem(nx,xn,'.');
ylabel('x(n)')subplot(2,2,2),stem(nh,hn,'.');ylabel('h(n)')
subplot(2,1,2);ny=1:L;stem(ny,real(yn),'.');ylabel('y(n)')
tic,yn=conv(xn,hn);toc

（4）用MATLAB语言编写程序，完成以下功能：读入一副图像，添加噪声，利用离散傅里叶变换频域分析方法，将图像在频域中进行噪声去除处理，将处理后的图像与原始图像相比

clear;
clc;
F = imread('D:\系统默认\图片\7efb984656f9000923eba0381a090a36_1.jpg');
F1 = rgb2gray(F);
subplot(2,2,1),imshow(F);
title('原色图像');
subplot(2,2,2),imshow(F1);
title('灰阶图像');
F1 = imnoise(F1,'gaussian',0.1);   %添加噪声
subplot(2,2,3),imshow(F1);
title('噪声图像');

FFT = fft2(F1);
myangle = angle(FFT);             %相位谱(没有进行移位的)
FS = abs(fftshift(FFT));          % 移位，使低频成分集中到图像中心，并得到幅度谱

[m,n] = size(FS);
FS(1:m/2-100,:) = 0;
FS(m/2+100:m,:) = 0;
FS(m/2-100:m/2+100,1:n/2-100) = 0;
FS(m/2-100:m/2+100,n/2+100:n) = 0;

aaa = ifftshift(FS);          % 将处理后的幅度图反移位，恢复到正常状态
bbb = aaa.*cos(myangle) + aaa.*sin(myangle).*1i;      % 幅度值和相位值重新进行结合，得到复数
fr = abs(ifft2(bbb));               % 进行傅里叶反变换，得到处理后的时域图像
ret = im2uint8(mat2gray(fr));       
subplot(2,2,4),imshow(ret);       %去除高频成分后的图像
title('去除噪声图像');

三、 实验中遇到的问题以及解决过程

在写代码时未注意中英文与半角输入法的区别导致代码无法正常运行

程序运行时，系统提示缺少某些函数库等问题。通过csdn下载相关插件解决了问题。

四、实验提问

（1）每行末尾……是什么意思？

代表这行内容还没有结束。

（2）= =是什么意思？

等于的意思。

五、实验总结

这次实验应用MATLAB数字信号处理工具箱处理离散信号（DFT），对各种典型序列做傅立叶变换，验证N点DFT的物理意义，严整频域采样与时域采样的对偶性及快速卷积算法。

通过本次内容让我对实验的基本内容了解的更加深刻，并发现了自己的基础错误，收获非常大。