应用FFT实现信号频谱分析(杭电)

信号、系统与信号处理实验Ⅱ

实验报告

实验名称：应用FFT实现信号频谱分析

一、实验目的

1、能够熟练掌握快速离散傅里叶变换FFT的原理及应用FFT进行频谱分析的基本方法。

2、对离散傅里叶变换的主要性质及FFT在数字信号处理中的重要作用有进一步的了解。

二、实验内容与要求

(2)编写一个调用FFT 函数的通用程序，可计算下列三种序列的离散频谱。

指数序列：v1(n) = 0.9^n*u(n);

周期为N 的余弦序列：v2(n) = cos(2*pi*n/N)，且0≤n≤N-1

复合函数序列：v3(n) = 0.9sin(2*pi*n/N) + 0.6sin(2*pi*n*3/N)

(3)计算实指数序列v1(n)的N点离散频谱V1(k),记录N为不同的2的幂次方时的V1(k)

值，并与理论值 V1(e^jwk )进行分析比较。

(4)计算周期为N 的余弦序列v2(n)的N 点FFT、2N 点FFT 及(N+2)点FFT,记录结果并

作为分析说明。

(5)已知信号x(t) = 0.15sin(2*pi*f1*t)+sin(2*pi*f2*t)-0.1sin(2*pi*f3*t),其中

f1 = 1Hz, f2 =2Hz,f3= 3Hz,取样频率为32Hz。编程实现：

① 32 点FFT，画出其幅度谱。

② 64 点FFT，画出其幅度谱，比较两者间的差异，思考实际频率与离散频谱图中

横坐标k的对应关系。

三、实验程序与结果

指数序列：v1(n) = 0.9^n*u(n);v1(n)的N点离散频谱V1(k),记录N为不同的2的幂次方时的V1(k)值，并与理论值 V1(e^jwk )进行分析比较。

clear all

N=16;

N1=16;

N2=32;

N3=64;

n=0:N-1;

xn=0.9.^n;

XK1=fft(xn,N1);

magXK1=abs(XK1);

phaXK1=angle(XK1);

XK2=fft(xn,N2); magXK2=abs(XK2); phaXK2=angle(XK2); XK3=fft(xn,N3); magXK3=abs(XK3); phaXK3=angle(XK3);

subplot(4,1,1);

plot(n,xn);

xlabel('n');

ylabel('x(n)');

title('x(n) ');

subplot(4,1,2);

k=0:length(magXK1)-1; stem(k,magXK1); xlabel('k');

ylabel('|X(k)|');

title('X(k) N1 = 16'); subplot(4,1,3);

k=0:length(magXK2)-1; stem(k,magXK2); xlabel('k');

ylabel('|X(k)|');

title('X(k) N2 = 32'); subplot(4,1,4);

k=0:length(magXK3)-1; stem(k,magXK3); xlabel('k');

ylabel('|X(k)|');

title('X(k) N3 = 64');

5

10

15

n

x (n )

k

|X (k )|

X (k) N1 = 16

5

10

15

20

25

30

35

k

|X (k )|

X (k) N2 = 32

10

20

30

40

50

60

70

k

|X (k )|

X (k) N3 = 64

周期为N 的余弦序列：v2(n) = cos(2*pi*n/N)，且0≤n ≤N-1；计算周期为N 的余弦序列v2(n)的N 点FFT 、2N 点FFT 及(N+2)点FFT,记录结果并作为分析说明。

clear all N=20; N1=20; N2=40; N3=22; n = 0:N-1;

xn =cos(2*pi*n/N);

XK1=fft(xn,N1); magXK1=abs(XK1); phaXK1=angle(XK1); XK2=fft(xn,N2); magXK2=abs(XK2); phaXK2=angle(XK2); XK3=fft(xn,N3); magXK3=abs(XK3); phaXK3=angle(XK3);

subplot(4,1,1);

plot(n,xn); xlabel('n'); ylabel('x(n)'); title('x(n) '); subplot(4,1,2);

k=0:length(magXK1)-1; stem(k,magXK1); xlabel('k');

ylabel('|X(k)|'); title('X(k) N = 20'); subplot(4,1,3);

k=0:length(magXK2)-1; stem(k,magXK2); xlabel('k');

ylabel('|X(k)|'); title('X(k) N = 40'); subplot(4,1,4);

k=0:length(magXK3)-1; stem(k,magXK3); xlabel('k');

ylabel('|X(k)|'); title('X(k) N = 22');

n

x (n )

k

|X (k )|

k

|X (k )|

X (k) N = 40

k

|X (k )|

X (k) N = 22

复合函数序列：v3(n) = 0.9sin(2*pi*n/N) + 0.6sin(2*pi*n*3/N)。 clear all N = 40; n = 0:N-1;

xn = 0.9*sin(2*pi.*n/N)+0.6*sin(2*pi*n*3/N); XK = fft(xn,N); magXK = abs(XK); phaXK = angle(XK); subplot(1,2,1); stem(n,xn); xlabel('n'); ylabel('x(n)'); title('x(n) N =40'); subplot(1,2,2);

k = 0:length(magXK)-1; stem(k,magXK); xlabel('k');

ylabel('|X(k)|'); title('X(k) N =40');

n

x (n )

x(n) N =40

k

|X (k )|

已知信号x(t) = 0.15sin(2*pi*f1*t)+sin(2*pi*f2*t)-0.1sin(2*pi*f3*t),其中f1=1Hz,f2=2Hz,f3=3Hz,取样频率为32Hz 。① 32 点FFT ，画出其幅度谱。② 64 点FFT ，画出其幅度谱，比较两者间的差异，思考实际频率与离散频谱图中横坐标k 的对应关系。

clear all N=32; N1=32; N2=64; n=0:N-1;

t=n/N;

xn=0.15*sin(2*pi*t)+sin(4*pi*t)-0.1*sin(6*pi*t); XK1=fft(xn,N1); magXK1=abs(XK1); phaXK1=angle(XK1); XK2=fft(xn,N2); magXK2=abs(XK2); phaXK2=angle(XK2); subplot(3,1,1); stem(n,xn); xlabel('n'); ylabel('x(n)'); title('x(n)'); subplot(3,1,2);

k=0:length(magXK1)-1; stem(k,magXK1); xlabel('k');

ylabel('|X(k)|'); title('X(k) 32'); subplot(3,1,3);

k=0:length(magXK2)-1; stem(k,magXK2); xlabel('k');

ylabel('|X(k)|'); title('X(k) 64');

051015

20253035

n x (n )

x(n)

k |X (k )|

X (k) 32

k

|X (k )|

X (k) 64

四、仿真结果分析

从上面32点FFT的频谱图中，可以看到6条谱线，其中后三条为前面对称的，即真正的分量有三个，对应原信号的三个频率。可以看出当频域取样点数增加后，可以使得原来看不见的频谱分量变得可以看得到。若对序列后面补零，也能看到之前看不到的频谱分量。可以通过增加抽样点数N，选择合适的窗函数来加以解决频谱混叠和频率泄漏。

五、实验问题解答与体会

1、利用DFT对连续信号进行傅里叶分析可能造成哪些误差以及造成这些误差的原因？

在运用DFT进行频谱分析的时候可能有三种误差：

（1）混淆现象：序列的频谱是采样信号频谱的周期延拓，当采样速率不满足Nyquist 定理，经过采样就会发生频谱混淆。

（2）泄露现象：实际中的信号序列往往很长，甚至是无线长序列。用截短的序列使用较短的DFT来对信号进行频谱分析。这种截短等价于给原始信号序列乘以一个矩

形窗函数，而矩形窗函数的频谱不是有限带宽的，从而它和原信号的频谱进行卷积以后会扩展原信号的频谱。

（3）栅栏效应：因为DFT是对单位圆上Z变换的均匀采样，所以它不可能将频谱视为一个连续函数。这样就产生了栅栏效应。

2、实验心得与体会

通过本次实验，我理解了快速离散傅里叶变换FFT的原理及应用FFT进行频谱分析的基本方法。而且对离散傅里叶变换的主要性质及FFT在数字信号处理中的重要作用有进一步的了解。