数字信号处理上机实验.doc

实验1 抽样定理的实验体会
实验内容：把下述三个连续时间信号转换成离散时间信号，在计算机上绘出的图形。为抽样频率。自行依次选取不同的抽样频率，如等。
（1） 工频信号：，，
Dt=0.00005;t=-0.005:Dt:0.05;
A=220;
fo=50;
xa=A*sin(2*pi*fo*t);
Ts=0.04;n=-25:1:25;
x=A*sin(2*pi*fo*n*Ts);
stem(n,x,'fill');
grid on;

图1.1 fs=25Hz时的图形

图1.2 fs=50Hz时的图形

图1.3 fs=100Hz时的图形

图1.3 fs=250Hz时的图形

（2） 衰减正弦信号：，，，
Dt=0.00005;t=-0.005:Dt:0.05;
A=2;a=0.5;fo=2;
xa=A*exp(-a*t).*sin(2*pi*fo*t);
Ts=1;n=-25:1:25;
x=A*exp(-a*n*Ts).*sin(2*pi*fo*n*Ts);
stem(n,x,'fill');
grid on;

图2.1 fs=1Hz时的图形

图2.2 fs=2Hz时的图形

图2.3 fs=4Hz时的图形

图2.4 fs=10Hz时的图形
（3） 谐波信号：，，，，
Dt=0.00005;t=-0.005:Dt:0.05;
A1=1;A2=0.5;A3=0.2;fo=5;
xa=A1*sin(2*pi*fo*t)+A2*sin(2*pi*fo*2*t)+A3*sin(2*pi*pi*3*t);
Ts=0.4;n=-25:1:25;
x=A1*sin(2*pi*fo*n*Ts)+A2*sin(2*pi*fo*2*n*Ts)+A3*sin(2*pi*pi*3*n*Ts);
stem(n,x,'fill');
grid on;

图3.1 fs=2.5Hz时的图
图3.2 fs=5Hz时的图形

图3.3 fs=10Hz时的图形

图3.4 fs=25Hz时的图形

实验2 离散信号的DTFT和DFT
实验内容： 分别计算16点序列 的16点和32点DFT，绘出幅度谱图形，并绘出该序列的DTFT图形。
实验要求：讨论DTFT和DFT之间的相互关系。说明实验产生的现象的原因。
N1=16;N2=32;
n1=0:N1-1;
n2=0:N2-1;
Xn=cos(5*pi*n/16);
x1k=fft(Xn,N1);
x2k=fft(Xn,N2);
subplot(2,1,1);stem(n1,abs(x1k),'.');axis([0,40,0,40]);ylabel('|x1(k)|')
title('16点的DFT[Xn]')
subplot(2,1,2);stem(n2,abs(x2k),'.');axis([0,40,0,40]);ylabel('|x2(k)|')
title('32点的DFT[Xn]')

图1 基本序列的离散傅里叶变换

num=[1 0.556 -0.383 -0.981 -0.707 0.195 0.924 0.831 0 -0.831 -0.924 -0.195 0.707 0.981 0.383 -0.556];
den=[1];
[h,w]=freqz(num,den,256,'whole',1);
plot(w,abs(h));

图2 序列的DTFT图形

DTFT和DFT之间的相互关系:DFT可以看做DTFT在区间［0, 2π］上的N点等间隔采样值, 其采样间隔为ωN=2π/N， 这就是DFT的物理意义。显而易见，DFT的变换区间长度N不同， 表示对X(ejω)在区间［0, 2π］上的采样间隔和采样点数不同， 所以DFT的变换结果也不同。

实验3 正弦信号抽样的实验
给定信号，现对x(t)抽样，设抽样点数N=16. 我们知道正弦信号的频谱是在处的函数，将x(t)抽样变成x(n)后，若抽样率及数据长度N取得合适，那么x(n)的DFT也应是在处的函数，由Parseval定理，有
表示x(n)的DFT在50Hz处的谱线，若上式不等，说明X(k)在频域有泄露。给定以下抽样频率（a），（a），（c），（1）分别得到x(n)及计算其X(k)，并用Parseval定理研究其泄露情况；
（2）当取，N=16时，在抽样点后面再补N个零，得到，这时是32点序列，求的DFT ，观察正弦信号补零的影响。
（3）观察抽样得到x(n)及X(k)，总结对正弦信号抽样应掌握的原则；
（1） f=1;
fs=2;
n1=[0:1:15];
f1fs/length(n1);
xa1=sin(2*pi*n1/fs);
for n=1:1:16
a(n)=xa1(n)*xa1(n);
t=a(n)
end
subplot(2,1,1)
stem(n1,xa1)
xlabel('n');ylabel('x(n)');
xk1=fft(xa1);
xk1=abs(xk1);
subplot(2,1,2)
stem(n1*f1,xk1);
xlabel'f');ylabel('X(k)');
当fs=2时，t= 0 1.4998e-32 5.9990e-32 1.3498e-31 2.3996e-31 3.7494e-31 5.3991e-31 7.3488e-31 9.5985e-31 1.2148e-30 1.4998e-30 2.4008e-29 2.1597e-30 3.8442e-30 2.9395e-30 2.9049e-29
同理可得fs=3,fs=4时的t值。经计算知当fs=2或fs=3时，Et不等于Ef, 由Parseval定理知X(k)在频域有泄漏；当fs=4时，Et等于Ef，X(k)不存在频域的泄漏。
\
图1 fs=2时X（n）及X（k）的图形

图2 fs=3时X（n）及X（k）的图形

图3 fs=4时X（n）及X（k）的图形
（2） M=1;
N1=length(n1);
z=zeros(1,M*N1);
N=(M+1)*N1;
n3=[0:1:N-1];
f3=8/N;
x=[xa1,z];
figure(2);
subplot( 2,2,1);
stem(n1,xa1);
xlabel('n');ylabel('x(n)');
subplot(2,2,2);
stem(n1*f1,xk1);
xlabel('f');ylabel('X(k)');
subplot(2,2,3)
stem(n3,x);
xlabel('n');ylabel('x(n)');
subplot(2,2,4);
xk=fft(x,N);
xk=abs(xk);
stem(n3*f3,xk)
xlabel('f');ylabel('X(k)

有图可知，对正弦信号抽样补零后X(k)发生了频域的泄漏。
（3）对正弦信号抽样应掌握的原则：
抽样频率应该为信号频率的整数倍，至少应取三倍，最好为四倍；抽样点数应包含整周期，一个周期内最好为四个点，以使数据点数为2的整数次幂；否则，尽管数据长度加大，反而泄漏更大。

实验4 快速Fourier变换(FFT)及其应用

实验中用到的信号序列：
a) Gaussian序列

b) 衰减正弦序列

c) 三角波序列

d) 反三角波序列

(1)、观察高斯序列的时域和幅频特性，固定信号xa(n)中参数p=8，改变q的值，使q分别等于2，4，8，观察它们的时域和幅频特性，了解当q取不同值时，对信号序列的时域幅频特性的影响；固定q=8，改变p，使p分别等于8，13，14，观察参数p变化对信号序列的时域及幅频特性的影响，观察p等于多少时，会发生明显的泄漏现象，混叠是否也随之出现？记录实验中观察到的现象，绘出相应的时域序列和幅频特性曲线。
p=8;
q=2;
N=16;
n=0:1:15;
x=exp(-(n-p).*(n-p)/q);
subplot(2,1,1);
stem(n,x,'fill');
x1=fft(x,N);
subplot(2,1,2);
stem(n,abs(x1),'.');

图1.1 p=8,q=2时的时域序列和幅频特性

图1.2 p=8,q=4时的时域序列和幅频特性

图1.3 p=8,q=8时的时域序列和幅频特性

图1.4 p=13,q=8时的时域序列和幅频特性

图1.5 p=14,q=8时的时域序列和幅频特性
由以上图形可知，1)当固定信号xa(n)中参数p=8，改变q的值使其增大时，信号的时域序列值增大，且均在n=P=8时取得最大值并且对称；幅频序列值也增大。2）当固定q=8，改变p使其增大时，信号的时域序列值右移且在n=p时大，幅频序列值减小。3）当p=8时会发生明显的泄漏现象，混叠也随之出现。
(2)、观察衰减正弦序列xb(n)的时域和幅频特性，a=0.1，f=0.0625，检查谱峰出现位置是否正确，注意频谱的形状，绘出幅频特性曲线，改变f，使f分别等于0.4375和0.5625，观察这两种情况下，频谱的形状和谱峰出现位置，有无混叠和泄漏现象？说明产生现象的原因。
N=16;
n=0:1:15;
a=0.1;
f=0.0625；
x=exp(-a*n).*sin(2*pi*f*n);
subplot(3,1,1);
stem(n,x,'fill');
x1=fft(x,N);
subplot(3,1,2);
stem(n,x1,'fill');
subplot(3,1,3);
stem(n,abs(x1),'.');

图2.1 a=0.1,f=0.0625的时域序列、DFT频谱、幅频特性

图2.2 a=0.1,f=0.4375时的时域序列、DFT频谱、幅频特性

图2.3 a=0.1,f=0.5625时的时域序列、DFT频谱、幅频特性
当a=0.1，f=0.0625时，谱峰出现位置正确。改变f，使f分别等于0.4375和0.5625时，由图可知，观察这两种情况下，频谱的形状和谱峰出现位置，出现了混叠和泄漏现象。
实际中我们往往用截短的序列来近似很长的甚至是无限长的序列，这样可以使用较短的DFT来对信号进行频谱分析，所得的频谱是原序列频谱的扩展。这个过程中产生了泄漏。而泄漏不能与混叠完全分开，因为泄漏导致频谱的扩展，从而造成混叠。为了减少泄漏的影响，可以选择适当的窗函数使频谱的扩散减至最小。避免混叠现象的唯一方法是保证采样速率足够高，使频谱混叠现象不致出现。
(3)、观察三角波和反三角波序列的时域和幅频特性，用N=8点FFT分析信号序列xc(n)和xd(n)的幅频特性，观察两者的序列形状和频谱曲线有什么异同？绘出两序列及其幅频特性曲线。
在xc(n)和xd(n)末尾补零，用N=16点FFT分析这两个信号的幅频特性，观察幅频特性发生了什么变化？两情况的FFT频谱还有相同之处吗？这些变化说明了什么？
N=8;
n = 0:7;
x = zeros(size(n))
mask = (n >= 0) & (n <= 3);
x(mask) = n(mask) + 1;
mask = (n >=4);
x(mask) = 8 - n(mask);
subplot(2,2,1);
stem(n,x,'fill');
subplot(2,2,3);
x1=fft(x,N);
stem(n,abs(x1),'.');
n1=0:7;
y=zeros(size(n));
mask1 = (n >= 0) & (n <= 3);
y(mask1) =4- n(mask1);
mask1 = (n >=4);
y(mask) = n(mask)-3;
subplot(2,2,2);
stem(n,y,'fill');
subplot(2,2,4);
y1=fft(y,N);
stem(n,abs(y1),'.');

N=8时，三角波和反三角波的序列形状不同，频谱曲线相同。
N=16时，频谱不再相同，说明了序列x(n)的长度N对DFT有影响。

图3.1 N=8时三角波和反三角波序列的时域和幅频特性
N=16;
n = 0:15;
x = zeros(size(n));
mask = (n >= 0) & (n <= 3);
x(mask) = n(mask) + 1;
mask = (n >=4)&(n<=7);
x(mask) = 8 - n(mask);
mask = (n >=8);
mask=0;
subplot(2,2,1);
stem(n,x,'fill');
subplot(2,2,3);
x1=fft(x,N);
stem(n,abs(x1),'.');
n1=0:15;
y=zeros(size(n));
mask1 = (n >= 0) & (n <= 3);
y(mask1) =4- n(mask1);
mask1 = (n >=4)&(n<=7);
y(mask1) = n(mask1)-3;
mask1 = (n >=8);

mask1=0;
subplot(2,2,2);
stem(n,y,'fill');
subplot(2,2,4);
y1=fft(y,N);
stem(n,abs(y1),'.');

图3.2 N=16时三角波和反三角波序列的时域和幅频特性
(4)、一个连续信号含两个频率分量，经采样得
x(n)=sin2π*0.125n+cos2π*(0.125+Δf)n n=0,1……,N-1
已知N=16，Δf分别为1/16和1/64，观察其频谱；当N=128时，Δf不变，其结果有何不同，为什么？
N=16;
f=1/16;
n=0:1:15;
x=sin(2*pi0.125*n)+cos(2*pi*(0.125+f)*n);
x1=fft(x,N);
stem(n,abs(x1),'.');

图4.1 N=16，Δf=1/16时的频谱

图4.2 N=16，Δf=1/64时的频谱

图4.3 N=128，Δf=1/16时的频谱

图4.4 N=128，Δf=1/64时的频谱
(5)、用FFT分别实现xa(n)（p＝8，q＝2）和 xb(n)（a＝0.1，f＝0.0625）的16点圆周卷积和线性卷积。
当N=16时，圆周卷积不等于线性卷积。要用FFT计算圆周卷积需在序列后面补零使其长度L=16+16-1=31。
M=16;
N=16;
nx=0:15;
nh=0:15;
xn=exp(-(nx-8).*(nx-8)/2);
hn=exp(-0.1*nh).*sin(2*pi*0.0625*nh);
L=pow2(nextpow2(M+N-1));
XK=fft(xn,L);
HK=fft(hn,L);
YK=XK.*HK;
yn=ifft(YK,L);
ny=0:L-1;
subplot(3,1,1);stem(nx,xn,'.');title('x(n)');
subplot(3,1,2);stem(nh,hn,'.');title('h(n)');
subplot(3,1,3);stem(ny,abs(yn),'.');title('y(n)');

图5 xa(n)、xb(n)及其线性卷积波形