数字信号处理实验总结.doc
《数字信号处理实验总结.doc》由会员分享，可在线阅读，更多相关《数字信号处理实验总结.doc（33页珍藏版）》请在冰点文库上搜索。
实验一离散信号及运算
一、实验目的
1.掌握MATLAB语言的基本功能及实现方法；
2.掌握MATLAB中各种常用序列的表示和显示方法；
3.熟练运用MATLAB进行离散信号的各种运算。
二、实验原理
我们所接触的信号大多为连续信号，而计算机及其他设备处理的大多为数字信号。
为了便于处理，往往要对信号进行处理使之变成离散数字信号。
对信号进行时间上的量化（即采样）是对信号作数字化处理的第一个环节，要求理解采样的原理和采样的性质，知道采样前后信号的变化及对离散信号和系统的影响。
三、实验内容
1、用MATLAB实现下列序列，并画出图形：
①单位采样序列移位，；
提示：
实现单位采样序列：
，可通过以下语句实现：
x=zeros（1,N）;x
（1）=1;
n=0:
10;
x=[zeros（1,3）,1,zeros（1,7）];
stem（n,x）；
②单位阶跃序列移位，
提示：
实现单位阶跃序列：
，可通过以下语句实现：
x=ones（1,N）;
n=0:
10;x=[zeros（1,3）,1,ones（1,7）];stem（n,x）
③正弦序列，，其中A=2；f=10；=0.005；
A=2;f=10;Ts=0.005;n=0:
10;
x=A*sin（2*pi*f*n*Ts）;stem（n,x）
③指数序列，
n=0:
10;x=0.9.^n;stem（n,x）
④复指数序列，，画出该序列的实部、虚部，幅值和相位。
提示：
可通过下列语句实现：
实部real（x），虚部imag（x），幅值abs（x），相位angle（x）
n=-20:
20;x=exp（0.05+j*pi/6*n）;
xr=real（x）;xi=imag（x）;xm=abs（x）;xa=angle（x）;
figure;
subplot（411）;stem（n,xr）;title（'实部'）;
subplot（412）;stem（n,xi）;title（'虚部'）;
subplot（413）;stem（n,xm）;title（'模'）;
subplot（414）;stem（n,xa）;title（'相角'）;
2、用MATLAB实现两个序列相加：
序列1：
x1=[10.50.30]，n1=1:
4；
序列2：
x2=[0.20.30.40.50.81]，n2=1:
6；实现x=x1+x2，n=1:
8，并画出x的图形。
提示：
MATLAB中可用算术运算符“+”实现序列相加，但两个序列的长度必须相等。
如果序列长度不等，或者长度虽然相等但采样的位置不同，就不能运用“+”了。
当两序列的长度不等或位置不对应时，首先应使两者位置对齐，然后通过zeros函数左右补零使其长度相等后再进行相加。
x1=[1,0.5,0.3,0];x2=[0.2,0.3,0.4,0.5,0.8,1];
n1=1:
4;n2=1:
6;n=1:
8;x3=[x1,0,0];x=x3+x2;
y=[x,0,0];stem（n,y）
3、用MATLAB实现序列的反转：
实现，序列x（n）采用，并画出y（n）的图形。
提示：
可利用fliplr（x）函数，例如：
x=[1234];y=fliplr（x）;结果为y=[4321]，要实现对fliplr（x）函数进行合理运用。
n=0:
10;
x=[zeros（1,3）,1,ones（1,7）];y=fliplr（x）;
n1=-fliplr（n）;stem（n1,y）;
4、序列的尺度变换，实现插值和抽取：
已知序列用MATLAB分别实现下列尺度变换。
提示：
可对序列x（n）的下标进行取余计算，余数为零即为插值和抹去的点，函数如下：
mod（nx,m），nx为序列x（n）的下标，m为插值或抽取的倍数。
clearall;
x=[1,2,3,4,5,6,7,8];
n=-4:
3;n1=n（1:
2:
length（n））;
y1=x（1:
2:
length（x））;subplot（211）;stem（n1,y1）;
axis（[-5307]）;y2=zeros（1,16）;
fork=1:
8
y2（2*k）=x（k）;
End;
subplot（2,1,2）;stem（-8:
6,y2（2:
end））;
四、思考题
1.若用C语言实现有限长序列的加法，编程如何实现？
与MATLAB相比，优缺点，繁简度如何？
2.MATLAB中可用算术运算符“+”实现序列相加，用算术运算符“*”实现序列相乘，试问用MATLAB求任意序列相加、相乘时，分别应注意什么？
3.分析stem（n,x）和plot（t,x）函数使用时的异同点。
实验二离散时间系统的表达及计算
一、实验目的
1.掌握离散系统的性质、输入输出关系及卷积运算；
2.掌握离散系统常系数线性差分方程的解法。
二、实验原理
一个离散时间系统在数学上的定义是将输入序列x（n）映射成输出序列y（n）的唯一性变
换或运算。
它的输入是一个序列，输出也是一个序列，其本质是将输入序列转变成输出序列的一个运算。
在数字信号处理中，通常研究的离散系统大都为LSI系统。
三、实验内容
1.求两个序列的卷积（离散线性卷积）
对于一个LSI系统，设其输入序列为，其中。
系统为单位冲激响应，编程求出输出序列，并画出图形。
提示：
输出序列应是输入序列和系统单位冲激相应序列的卷积，在Matlab中可用conv函数来实现两个序列的卷积。
但是使用conv函数求卷积时，两序列卷积后的输出序列的下标不是任意的（具体是多少，查函数），所以需要自己构造下标。
n=-5:
50;
u1=stepseq（0,-5,50）;u2=stepseq（10,-5,50）;
x=u1-u2;h=（（0.9）.^n）.*u1;
[y,ny]=conv1（x,n,h,n）;stem（ny,y）;
function[x,n]=stepseq（n0,n1,n2）;
ifnargin~=3
disp（'Usage:
Y=stepseq（n0,n1,n2）'）;
elseif（（n0n2）|（n1>n2））
error（'arguments?
must?
satisfy?
n1<=n0<=n2'）;
End;
n=[n1:
n2];x=[（n-n0）>=0];
function[y,ny]=conv1（x,nx,h,nh）
ifnargin~=4
disp（'Usage:
Y=conv_m（x,nx,h,nh）'）;
return;end;
nyb=nx
（1）+nh
（1）;nye=nx（length（x））+nh（length（x））;
ny=[nyb:
nye];y=conv（x,h）;
2.求两个序列的卷积（离散周期卷积）
已知序列，，计算它们的周期卷积。
提示：
将和的主值序列x（n）和y（n）作线性卷积，然后再将线性卷积序列以N为周期进行周期延拓。
clearall;
x=[1,0,1,2,1];y=[1,1,0,1,2];
z=conv（x,y）;
m=0:
4;
h=[z（6:
9）,0];h1=h+z（1:
5）;h2=[h1h1h1];
stem（h2）
3.已知LSI系统的差分方程为：
，其中
求单位冲激响应，并画出图形。
提示：
可以使用filter函数来完成，其调用格式为，注意b，a，x参数都代表什么，如何设置。
clearall
b=[1];a=[1-10.9];
n=-20:
100;N=120;
x=[zeros（1,20）,1,zeros（1,100）];y=filter（b,a,x）;
stem（n,y）;axis（[-20100-22]）;
4.已知LSI系统的差分方程为：
，其中
求单位阶跃响应，并画出图形。
提示：
可以使用filter函数来完成，其调用格式为，注意b，a，x参数都代表什么，如何设置。
b=[1];a=[1-10.9];
N=120;n=-20:
100;
x=[zeros（1,20）,1,ones（1,100）];
g=filter（b,a,x）;stem（n,g）;
思考题
1.conv函数来实现两个序列的卷积，其默认下标为？
（10分）
2.离散线性卷积与离散周期卷积的区别？
（20分）
3.解差分方程时，a、b分别代表什么？
（10分）
当差分方程为时，a、b分别为多少？
4.解释filter函数？
（10分）
实验三z变换及反变换
一．实验目的
1.通过Matlab编程，熟悉z变换定义及逆z变换常用方法，加深对z变换性质及其收敛域
的理解；
2.通过Matlab编程实现离散信号及系统的z域分析；
3.掌握利用Matlab编程求解差分方程的方法。
二．实验原理
z变换性质、求逆z变换的方法、系统的z域表示方法及差分方程的求解方法。
三．实验内容
1.z正变换：
已知序列x1=[1,0,5,6]，其中；与x2=[1,3,5,2]，其中，，求其卷积信号x=x1*x2的z变换。
提示：
例将序列x1进行Z变换，则根据公式
（与序列x1比较）
也就是说，没必要在Matlab命令窗口中显示出完整序列z变换的表达式，通过序列的位置信息和该位置信息上所对应的序列值，即可在报告中将该序列的z变换写出来。
clearall
x1=[1,0,5,6];x2=[1,3,5,2];
symsz;x=conv（x1,x2）;s=0;
forn=1:
length（x）
s=s+x（n）*z^（-n+1）;
End;s
x=131023434012
s=3/z+10/z^2+23/z^3+43/z^4+40/z^5+12/z^6+1
2.z反变换：
已知与，求的逆z变换；
提示：
思路与上题相同。
通过X1（z）和X2（z）找出序列x1（n）和x2（n），直接求x1（n）和x2（n）的卷积，即是求X1（z）·X2（z）的逆z变换。
注意序列x1（n）和x2（n）的n的取值范围。
n1=-1:
1;n2=-2:
1;
x1=[1,2,3];x2=[2,4,3,5];
y=conv（x1,x2）;n=n1
（1）+n2
（1）:
n1（3）+n2（4）;
stem（n,y）;
3.离散时间线性时不变系统的时域实现：
已知因果LSI系统的差分方程为：
，写出系统传递函数H（z）及其收敛域，编程实现以下内容：
①系统冲激响应h（n），并绘图；②系统的单位阶跃响应g（n），并绘图；
③系统函数的幅频响应H（ejw），并绘出幅频和相频特性。
提示：
①冲激响应可用impz（b,a,n）函数，可查看MatlabHelp了解它的相关用法。
impz函数：
数字滤波器的冲激响应，用来求系统的冲激响应h（n）。
b是差分方程中输入的系数，a是差分方程中输出的系数。
注意缺项的系数前要添加系数0。
这个函数本身就有画图功能。
如果n是一个整形变量，impz计算的是在这些整数点位置上的冲激响应，且从0开始。
②单位阶跃响应可用stepz（b,a,n）函数，可查看MatlabHelp了解它的相关用法。
与impz用法相似。
③幅频响应可用freqz（b,a,w），其中w=[0:
1:
500]*pi/500和freqz（b,a,n）分别求，可查看MatlabHelp了解它的相关用法。
①：
b=[10-1];a=[10-0.81];impz（b,a,100）;
②：
b=[10-1];a=[10-0.81];stepz（b,a,100）;
③：
a=[1,0,-0.81];b=[1,0,-1];freqz（b,a,50,'whole'）;
a=[1,0,-0.81];b=[1,0,-1];freqz（b,a,[0:
1:
500]*pi/500）;
四．思考题
1.画图语句后加grid，起什么作用？
（10分）
2.画图时，不加axis语句，为什么也可以画出图形？
此时横纵轴坐标怎样？
（10分）
3.freqz（b,a,w）中，w=[0:
1:
500]*pi/500，若改为：
freqz（b,a,n），n=[0:
1:
50]，会出现什么情况？
为什么？
（10分）
4.根据实验结果，写出系统函数的表达式H（z），并写出收敛域。
（10分）
5.系统稳定性和因果性的判决条件？
程序实现的主要思路？
（10分）
实验四离散时间傅立叶变换DFT（离散时间信号的频域分析）
五、实验目的
4.深刻理解有限长序列的离散傅立叶变换的概念、计算及意义；
5.掌握从离散傅立叶级数到离散傅立叶变换的过程；
6.掌握离散傅立叶变换的性质，熟悉信号的频域转换。
六、实验原理
一个信号的时域转换可以通过傅立叶变换（DFT）来完成。
有限长序列的DFT是其Z
变换在单位圆上的等距采样，或者说是序列傅立叶变换的等距采样，因此可以用于序列的频谱分析。
为了更好地理解DFT的概念，我们先了解周期序列及其离散傅立叶级数（DFS）表示。
七、实验内容
1.若是一个N=12的有限长序列，利用MATLAB计算它的DFT并画出图形。
提示：
在计算机上实现信号的频谱分析及其他方面的处理工作时，对信号的要求是：
在时域和频域都应该是离散的，且都应是有限长的，只有DFS（离散周期序列的傅立叶级数）在时域和频域都是离散的。
所以，DFT实际上来是自于DFS的，只不过仅在时域、频域各取一个周期而已。
本题要求自己编写一个函数，其中xn代表序列，N代表有限长序列的点数。
函数的算法采用序列求DFT变换的基本公式，即：
，其中，。
由于X（k）为复数，所以画图时，应画X（k）的幅度特性abs（）和相位特性angle（）。
clc
N=12;n=0:
N-1;k=0:
N-1;
x=cos（n*pi/6）;w=exp（-j*2*pi/N）;
A=n'*k;X=x*w.^A;X1=abs（X）;X2=angle（X）;
subplot（121）;stem（n,X1）;gridon;
subplot（122）;stem（n,X2）;gridon;
2.对于实序列
①分解成偶部和奇部；
②验证实序列的性质：
提示：
任意一个实序列都可以表示成偶对称序列与奇对称序列之和，即：
。
令，，这样分别求和，再求，将得到的X（k）取其和，与上述结果相比较，验证实序列的上述性质。
本题要求自己编写一个circevod（x）函数，它的功能是将一个实序列分解成与之和，并求出与。
其中需要有判断语句，判断这个序列是否是实序列；的求法可参考mod（-n,N）函数。
clearall;
n=0:
19;N=length（n）;m=mod（-n,N）+1;
x1=0.9.^n;x2=x1（m）;xe=1/2*（x1+x2）;xo=1/2*（x1-x2）;
subplot（321）;stem（n,xe）;title（'x1（n）的偶部'）;
subplot（322）;stem（n,xo）;title（'x1（n）的奇部'）;
Xe=dft（xe,N）;Xo=dft（xo,N）;
X=dft（x1,N）;ReX=real（X）;ImX=imag（X）;
subplot（323）;stem（n,Xe）;title（'xe（n）的DFT'）;
subplot（324）;stem（n,ReX）;title（'X的实部'）;
subplot（325）;stem（n,imag（Xo））;title（'xo的DFT'）;
subplot（326）;stem（n,ImX）;title（'X的虚部'）;
离散傅里叶变换的应用
1用DFT计算线性卷积：
（了解）
算法流程如下：
数字信号处理的优势是“实时实现”，即信号进来经处理后马上输出去。
然而：
引入两个问题：
lx（n）没有全部进入，如何实现卷积？
l全部进入再卷积，又如何保证实时实现？
解决方法通常是将较长序列进行分段，然后计算每段子序列和较短序列的线性卷积，最后再将各段线性卷积结果序列进行相加得到结果。
这类方法包括：
重叠相加法和重叠保留法。
2用DFT对模拟信号进行谱分析
编程实现：
给定模拟信号，以毫秒进行取样，
用DFT进行信号频谱分析：
①若要能区别该信号中的两个频率分量，试问取样信号的长度至少为多少？
共取多少采样点？
②所用DFT的点数N分别等于128、256、1024，画出信号的N点DFT的幅度谱。
③讨论幅度谱结果，N为多少时能分辨出信号中的所有频率分量？
clear;
N1=64;n1=0:
N1-1;N2=256;n2=0:
N2-1;N3=1024;n3=0:
N3-1;
t=0.01*n1;
x=2*sin（4*pi*t）+5*cos（8*pi*t）;X=dft（x,N1）;
subplot（321）;plot（t,x）;title（'x（t）---N=64'）;
subplot（322）;plot（t,abs（X））;title（'abs（X）---N=64'）;
t=0.01*n2;
x=2*sin（4*pi*t）+5*cos（8*pi*t）;X=dft（x,N2）;
subplot（323）;plot（t,x）;title（'x（t）---N=256'）;
subplot（324）;plot（t,abs（X））;title（'abs（X）---N=256'）;
t=0.01*n3;
x=2*sin（4*pi*t）+5*cos（8*pi*t）;X=dft（x,N3）;
subplot（325）;plot（t,x）;title（'x（t）---N=1024'）;
subplot（326）;plot（t,abs（X））;title（'abs（X）---N=1024'）;
八、思考题
1.圆周卷积与线性卷积的关系，二者是否相同？
在什么条件下相同？
2.离散信号的卷积运算有何用途？
3.DFT变换中补零对信号频谱的分辨率有影响吗？
补零后频谱图有何变化？
4.实验内容3中，若要能区别该信号中的两个频率分量，取样信号的长度至少为多少？
N分别为60、500时能否分辨出信号中所有频率分量？
与
（2）的结果比较。
实验五离散傅立叶变换和快速FFT
一．实验目的
1.在理论学习的基础上，通过本实验，进一步加深对DFT的算法原理及性质的理解（因为FFT只是DFT的一种快速算法，所以FFT的结果必然满足DFT的基本性质）；
2.熟悉并掌握按时间/频率抽取FFT算法原理和子程序的应用；
3.学习用FFT对连续信号和时域离散信号进行谱分析的方法；
4.了解应用FFT进行信号频谱分析过程中可能出现的问题（如：
混叠、泄露、栅栏效应等），
分析其原因，以便在实际中正确应用FFT。
二．实验原理
1.FFT原理：
DFT如果直接计算的话，计算量非常大，而且不利于信号的实时处理，在实际应用中遇到很大的困难，由此出现了很多快速的计算DFT的方法，在此我们以基2的时间抽取快速傅立叶算法为例。
2.混叠：
序列的频谱是被采样信号的周期延拓，当采样速率不满足Nyquist定理时，就会发生频谱混叠，使得采样后的信号序列频谱不能真实的反映原信号的频谱，避免混叠现象的唯一方法是保证采样速率足够高。
3.泄露：
用截短的序列来近似很长甚至是无限长的序列，这样可以使用较短的DFT对信号进行频谱分析，所得的频谱是原序列频谱的扩展。
为了减少泄露的影响，可以选择适当的窗函数使频谱的扩散减至最小。
4.栅栏效应：
DFT是对单位圆上Z变化的均匀采样，所以它不可能将频谱视为一个连续函数，用DFT来观察频谱就好像通过一个栅栏来观看一个图景一样，只能在离散点上看到真实的频谱，这样一些频谱的峰点或谷点被“栅栏”所拦住，不能被观察到。
减小栅栏效应的方法就是借助于在原序列的末端填补一些零值，从而变动DFT的点数，这实际上是人为地改变了对真实频谱采样的点数和位置，相当于搬动了每一根“栅栏”的位置，使频谱的峰点或谷点暴露出来。
三．实验内容
1.求序列x=[52741139]的快速傅立叶变换，并绘出幅频特性曲线；
提示：
采用fft（x,N）函数，N为FFT的点数，当序列长度小于N时，系统自动在序列末尾补零；当序列长度大于N时，系统自动截断序列多余的部分。
x=[52741139];
N=length（x）;xk=fft（x,N）;n=0:
N-1;
stem（n,abs（xk））;
2.设一序列中含有两种频率成分，f=2Hz，f=2.05Hz，采样率取为fs=10Hz，即
x（n）=sin（2fn/f）+sin（2fn/f），根据公式<，要区分出这两种频率成分，N必须满足多少？
⑴取x（n）（0n<128）时，计算128点FFT，并绘出幅频特性曲线；
⑵取x（n）（0n<128）时，补384个0，计算512点FFT，并绘出幅频特性曲线；
⑶取x（n）（0n<512）时，计算512点FFT，并绘出幅频特性曲线；
⑷分别改变FFT变换的点数和采样时间，对结果进行分析，对比看有何不同。
提示：
补零时可采用zeros函数；（4）中，N，t的取值为N=512，t=0.1n；N=1024，t=0.1n；N=512，t=0.05n；N=1024，t=0.05n四种情况。
f1=2;f2=2.05;fs=10;
N1=128;N2=512;
n1=0:
N1-1;n2=0:
N2-1;
n3=n1*fs/N1;n4=n2*fs/N2;n5=n2*fs/N2;
xn1=sin（2*pi*f1.*n1/fs）+sin（2*pi*f2.*n1/fs）;
xn2=sin（2*pi*f1.*n2/fs）+sin（2*pi*f2.*n2/fs）;
xk1=fft（xn1,N1）;xk2=fft（xn1,N2）;xk3=fft（xn2,N2）;
subplot（311）;plot（n3,abs（xk1））;title（'128点'）;
subplot（312）;plot（n4,abs（xk2））;title（'512点补零'）;
subplot（313）;plot（n5,abs（xk3））;title（'512点'）;
f1=2;f2=2.05;fs1=10;fs2=20;
升级会员 