信号处理实验指导书.docx
《信号处理实验指导书.docx》由会员分享,可在线阅读,更多相关《信号处理实验指导书.docx(53页珍藏版)》请在冰豆网上搜索。
信号处理实验指导书
《信号处理基础》实验指导书
计算机科学与技术学院
淮北师范大学2012
I
目录
实验
实验一线性时不变系统的响应............................................................................................1
实验二连续时间信号的频域分析.....................................................................................................6
实验三离散信号与系统的Z变换分析............................................................................................10
附录
附录一MATLAB环境.....................................................................................................................14
附录三基本绘图命令.....................................................................................................................18
附录四多项式的求值、求根和部分分式展开.............................................................................21
附录五符号积分变换.....................................................................................................................22
附录六信号与系统分析常用函数.................................................................................................24
信号处理基础实验指导书——实验部分
1
实验一线性时不变系统的响应
一、实验目的
1.熟悉连续时间系统的单位冲激响应、阶跃响应的意义及求解方法
2.熟悉连续(离散)时间系统在任意信号激励下响应的求解方法
3.熟悉应用MATLAB实现求解系统响应的方法
二、实验原理及实验内容
1.连续时间系统
对于连续的线性时不变系统(LTI)系统,当系统输入为f(t),输出为y(t),则输入与输出之间满足如
下的线性常系数微分方程:
()()
00
()()
nm
ij
ij
ij
aytbft
==
Σ=Σ,当系统输入为单位冲激信号δ(t)时产生的零状态响
应称为系统的单位冲激响应,用h(t)表示。
若输入为单位阶跃信号ε(t)时,系统产生的零状态响应则称为系
统的单位阶跃响应,记为g(t),如下图所示。
系统的单位冲激响应h(t)包含了系统的固有特性,它是由系统本身的结构及参数所决定的,与系统的
输入无关。
我们只要知道了系统的冲激响应,即可求得系统在不同激励下产生的响应。
因此,求解系统的
冲激响应h(t)对我们进行连续系统的分析具有非常重要的意义。
在MATLAB中有专门用于求解连续系统冲激响应和阶跃响应,并绘制其时域波形的函数impulse()
和step()。
如果系统输入为f(t),冲激响应为h(t),系统的零状态响应为y(t),则有:
y(t)=h(t)∗f(t)。
若已知系统的输入信号及初始状态,我们便可以用微分方程的经典时域求解方法,求出系统的响应。
但是对于高阶系统,手工计算这一问题的过程非常困难和繁琐。
在MATLAB中,应用lsim()函数很容易就能对上述微分方程所描述的系统的响应进行仿真,求出系
统在任意激励信号作用下的响应。
lsim()函数不仅能够求出连续系统在指定的任意时间范围内系统响应的
数值解,而且还能同时绘制出系统响应的时域波形图。
以上各函数的调用格式如下:
⑴impulse()函数
函数impulse()将绘制出由向量a和b所表示的连续系统在指定时间范围内的单位冲激响应h(t)的时域波
形图,并能求出指定时间范围内冲激响应的数值解。
impulse(b,a)以默认方式绘出由向量a和b所定义的连续系统的冲激响应的时域波形。
impulse(b,a,t0)绘出由向量a和b所定义的连续系统在0~t0时间范围内冲激响应的时域波形。
impulse(b,a,t1:
p:
t2)绘出由向量a和b所定义的连续系统在t1~t2时间范围内,并且以时间间隔p均
匀取样的冲激响应的时域波形。
y=impulse(b,a,t1:
p:
t2)只求出由向量a和b所定义的连续系统在t1~t2时间范围内,并且以时间间
隔p均匀取样的冲激响应的数值解,但不绘出其相应波形。
信号处理基础实验指导书——实验部分
2
⑵step()函数
函数step()将绘制出由向量a和b所表示的连续系统的阶跃响应,在指定的时间范围内的波形图,并
且求出数值解。
和impulse()函数一样,step()也有如下四种调用格式:
step(b,a)
step(b,a,t0)
step(b,a,t1:
p:
t2)
y=step(b,a,t1:
p:
t2)
上述调用格式的功能和impulse()函数完全相同,所不同只是所绘制(求解)的是系统的阶跃响应g(t),
而不是冲激响应h(t)。
⑶lsim()函数
根据系统有无初始状态,lsim()函数有如下两种调用格式:
①系统无初态时,调用lsim()函数可求出系统的零状态响应,其格式如下:
lsim(b,a,x,t)绘出由向量a和b所定义的连续系统在输入为x和t所定义的信号时,系统零状态响应
的时域仿真波形,且时间范围与输入信号相同。
其中x和t是表示输入信号的行向量,t为表示输入信号时
间范围的向量,x则是输入信号对应于向量t所定义的时间点上的取样值。
y=lsim(b,a,x,t)与前面的impulse和step函数类似,该调用格式并不绘制出系统的零状态响应曲线,
而只是求出与向量t定义的时间范围相一致的系统零状态响应的数值解。
②系统有初始状态时,调用lsim()函数可求出系统的全响应,格式如下:
lsim(A,B,C,D,e,t,X0)绘出由系数矩阵A,B,C,D所定义的连续时间系统在输入为e和t所定义的信号
时,系统输出函数的全响应的时域仿真波形。
t为表示输入信号时间范围的向量,e则是输入信号e(t)对应
于向量t所定义的时间点上的取样值,X0表示系统状态变量X=[x1,x2,…..xn]'在t=0时刻的初值。
[Y,X]=lsim(A,B,C,D,e,t,X0)不绘出全响应波形,而只是求出与向量t定义的时间范围相一致的系统
输出向量Y的全响应以及状态变量X的数值解。
显然,函数lsim()对系统响应进行仿真的效果取决于向量t的时间间隔的密集程度,t的取样时间间隔
越小则响应曲线越光滑,仿真效果也越好。
说明:
(1)当系统有初始状态时,若使用lsim()函数求系统的全响应,就要使用系统的状态空间描述法,即首先
要根据系统给定的方式,写出描述系统的状态方程和输出方程。
假如系统原来给定的是微分方程或系统函
数,则可用相变量法或对角线变量等方法写出系统的状态方程和输出方程。
其转换原理如前面实验四所述。
(2)显然利用lsim()函数不仅可以分析单输入单输出系统,还可以分析复杂的多输入多输出系统。
1-1:
若某连续系统的输入为e(t),输出为r(t),系统的微分方程为
y''(t)+5y'(t)+6y(t)=3f'(t)+2f(t)
①求该系统的单位冲激响应h(t)及其单位阶跃响应g(t)。
②若f(t)=e−2tε(t)求出系统的零状态响应y(t)
分析:
①求冲激响应及阶跃响应的MATLAB程序:
a=[156];b=[32];
subplot(2,1,1),impulse(b,a,4)
subplot(2,1,2),step(b,a,4)
画出运行结果图。
信号处理基础实验指导书——实验部分
3
②求零状态响应的MATLAB程序:
a=[156];b=[32];
p1=0.01;%定义取样时间间隔为0.01
t1=0:
p1:
5;%定义时间范围
x1=exp(-2*t1);%定义输入信号
lsim(b,a,x1,t1),%对取样间隔为0.01时系统响应进行仿真
holdon;%保持图形窗口以便能在同一窗口中绘制多条曲线
p2=0.5;%定义取样间隔为0.5
t2=0:
p2:
5;%定义时间范围
x2=exp(-2*t2);%定义输入信号
lsim(b,a,x2,t2),holdoff%对取样间隔为0.5时系统响应进行仿真并解除保持
画出运行结果图
1-2已知一个过阻尼二阶系统的状态方程和输出方程分别为:
010
'()()()
232
xtXtft
⎡⎤⎡⎤
=⎢⎥+⎢⎥⎣−−⎦⎣⎦
,r(t)=[01]X(t)。
若系统初始状态为X(0)=[4-5]T,求系统在f(t)=3e−4tε(t)作用下的全响应。
求全响应程序如下:
A=[01;-2-3];B=[02]';C=[01];D=[0];
X0=[4-5]';%定义系统初始状态
t=0:
0.01:
10;
E=[3*exp(-4*t).*ones(size(t))]';%定义系统激励信号
[r,x]=lsim(A,B,C,D,E,t,X0);%求出系统全响应的数值解
plot(t,r)%绘制系统全响应波形
画出运行结果图。
2.离散时间系统
LTI离散系统中,其输入和输出的关系由差分方程描述:
00
()()
nm
ij
ij
aykibfkj
==
Σ+=Σ+(前向差分方程)
00
()()
nm
ij
ij
aykibfknj
==
Σ−=Σ−+(后向差分方程)
当系统的输入为单位序列δ(k)时产生的零状态响应称为系统的单位函数响应,用h(k)表示。
当输入为ε(k)
时产生的零状态响应称为系统的单位阶跃应,记为:
g(k),如下图所示。
如果系统输入为e(k),冲激响应为h(k),系统的零状态响应为y(k),则有:
y(k)=h(k)∗f(k)。
与连
续系统的单位冲激响应h(t)相类似,离散系统的单位函数响应h(k)也包含了系统的固有特性,与输入序列
信号处理基础实验指导书——实验部分
4
无关。
我们只要知道了系统的单位函数响应,即可求得系统在不同激励信号作用下产生的响应。
因此,求
解系统的单位函数响应h(k)对我们进行离散系统的分析也同样具有非常重要的意义。
MATLAB中为用户提供了专门用于求解离散系统单位函数响应,并绘制其时域波形的函数impz()。
同样也提供了求离散系统响应的专用函数filter(),该函数能求出由差分方程所描述的离散系统在指定时间
范围内的输入序列作用时,产生的响应序列的数值解。
当系统初值不为零时,可以使用dlsim()函数求出
离散系统的全响应,其调用方法与前面连续系统的lsim()函数相似。
另外,求解离散系统阶跃响应可以通
过如下两种方法实现:
一种是直接调用专用函数dstep(),其调用方法与求解连续系统阶跃响应的专用函数
step()的调用方法相似;另一种方法是利用求解离散系统零状态响应的专用函数filter(),只要将其中的激
励信号看成是单位阶跃信号ε(k)即可。
函数的调用格式分别如下:
⑴impz()函数
impz(b,a)以默认方式绘出由向量a和b所定义的离散系统单位函数响应的时域波形。
impz(b,a,n)绘出由向量a和b所定义的离散系统在0~n(n必须为整数)的离散时间范围内单
位函数响应的时域波形。
impz(b,a,n1:
n2)绘出由向量a和b所定义的离散系统在n1~n2(n1、n2必须为整数)的离散
时间范围内单位函数响应的时域波形。
y=impz(b,a,n1:
n2)求出由向量a和b所定义的离散系统在n1~n2(n1、n2必须为整数)的离
散时间范围内单位函数响应的数值解,但不绘出波形。
⑵filter()函数
filter(b,a,x)其中a和b与前面相同,x是包含输入序列非零样值点的的行向量。
此命令将求出
系统在与x的取样时间点相同的输出序列样值。
1-3:
已知描述离散系统的差分方程为:
y(k)−0.25y(k−1)+0.5y(k−2)=f(k)+f(k−1),且已知系统
输入序列为1
2f(k)=()kε(k),
①求出系统的单位函数响应h(k)在-3~10离散时间范围内响应波形。
②求出系统零状态响应在0~15区间上的样值;并画出输入序列的时域波形以及系统零状态响应的波形
分析:
①求系统的单位函数响应的MATLAB程序:
a=[1,-0.25,0.5];b=[1,1,0];
impz(b,a,-3:
10),title('单位响应')%绘出单位函数响应在-3~10区间上的波形
②求零状态响应的MATLAB程序:
a=[1,-0.25,0.5];b=[1,1,0]
k=0:
15;%定义输入序列取值范围
x=(1/2).^k;%定义输入序列表达式
y=filter(b,a,x)%求解零状态响应样值
subplot(2,1,1),stem(k,x)%绘制输入序列的波形
title('输入序列')
subplot(2,1,2),stem(k,y)%绘制零状态响应的波形
title('输出序列')
列出运行结果,并画出波形图。
信号处理基础实验指导书——实验部分
5
三、预习要求
1.熟悉系统响应的求解方法
2.了解MATLAB语言中关于系统分析的各个函数如:
impulse、step、lsim、impz、filter等函数的调用方法:
四、实验报告要求
1.理论上计算出系统的单位冲激响应/单位函数响应、阶跃响应、零状态响应、全响应的表达式,并写出解
题过程。
2.记录仿真结果(包括数据和波形)。
3.写出程序清单。
4.实验总结(收获及体会)
信号处理基础实验指导书——实验部分
6
实验二连续时间信号的频域分析
一、实验目的
1.熟悉傅里叶变换的性质
2.熟悉常见信号的傅里叶变换
3.了解傅里叶变换的MATLAB实现方法
二、实验原理及实验内容
傅里叶变换是信号分析的最重要的内容之一。
从已知信号f(t)求出相应的频谱函数F(jω)的数学表
示为:
F(jω)∞f(t)e−jωtdt
−∞
=∫
f(t)的傅里叶变换存在的充分条件是f(t)在无限区间内绝对可积,即f(t)满足下式:
f(t)dt∞
−∞
∫<∞
但上式并非傅里叶变换存在的必要条件。
在引入广义函数概念之后,使一些不满足绝对可积条件的函数也
能进行傅里叶变换。
傅里叶反变换的定义为:
()1()
2
ftFjωejωtdω
π
∞
−∞
=∫。
在这一部分的学习中,大家都体会到了这种数学运算的麻烦。
在MATLAB语言中有专门对信号进行正
反傅里叶变换的语句,使得傅里叶变换很容易在MATLAB中实现。
在MATLAB中实现傅里叶变换的方
法有两种,一种是利用MATLAB中的SymbolicMathToolbox提供的专用函数直接求解函数的傅里叶变
换和傅里叶反变换,另一种是傅里叶变换的数值计算实现法。
下面分别介绍这两种实现方法的原理。
1.直接调用专用函数法
①在MATLAB中实现傅里叶变换的函数为:
F=fourier(f)对f(t)进行傅里叶变换,其结果为F(w)
F=fourier(f,v)对f(t)进行傅里叶变换,其结果为F(v)
F=fourier(f,u,v)对f(u)进行傅里叶变换,其结果为F(v)
②傅里叶反变换
f=ifourier(F)对F(w)进行傅里叶反变换,其结果为f(x)
f=ifourier(F,U)对F(w)进行傅里叶反变换,其结果为f(u)
f=ifourier(F,v,u)对F(v)进行傅里叶反变换,其结果为f(u)
由于MATLAB中函数类型非常丰富,要想了解函数的意义和用法,可以用mhelp命令。
如在命令窗
口键入:
mhelpfourier回车,则会得到fourier的意义和用法。
注意:
信号处理基础实验指导书——实验部分
7
(1)在调用函数fourier()及ifourier()之前,要用syms命令对所有需要用到的变量(如t,u,v,w)等进行说
明,即要将这些变量说明成符号变量。
对fourier()中的f及ifourier()中的F也要用符号定义符sym将其说
明为符号表达式。
(2)采用fourier()及fourier()得到的返回函数,仍然为符号表达式。
在对其作图时要用ezplot()函数,而
不能用plot()函数。
(3)fourier()及fourier()函数的应用有很多局限性,如果在返回函数中含有δ(ω)等函数,则ezplot()函数
也无法作出图来。
另外,在用fourier()函数对某些信号进行变换时,其返回函数如果包含一些不能直接表
达的式子,则此时当然也就无法作图了。
这是fourier()函数的一个局限。
另一个局限是在很多场合,尽管
原时间信号f(t)是连续的,但却不能表示成符号表达式,此时只能应用下面介绍的数值计算法来进行傅氏
变换了,当然,大多数情况下,用数值计算法所求的频谱函数只是一种近似值。
2-1求门函数f(t)=ε(t+1)−ε(t−1)的傅里叶变换,并画出幅度频谱图
MATLAB程序如下:
symstw%定义两个符号变量t,w
Gt=sym('Heaviside(t+1)-Heaviside(t-1)');%产生门宽为2的门函数
Fw=fourier(Gt,t,w);%对门函数作傅氏变换求F(jw)
FFw=maple('convert',Fw,'piecewise');%数据类型转换,转为分段函数,此处可以去掉
FFP=abs(FFw);%求振幅频谱|F(jw)|
ezplot(FFP,[-10*pi10*pi]);grid;%绘制函数图形,并加网格
axis([-10*pi10*pi02.2])%限定坐标轴范围
写出运行结果并画出波形图。
运行结果参考:
Fw=exp(i*w)*(pi*Dirac(w)-i/w)-exp(-i*w)*(pi*Dirac(w)-i/w)
%Dirac(w)为δ(ω),即傅立叶变换结果中含有奇异函数,故绘图前需作函数类型转换
FFw=-i*exp(i*w)/w+i*exp(-i*w)/w%FFw为复数
FFP=abs(-i*exp(i*w)/w+i*exp(-i*w)/w)%求FFw的模值
2-2求函数
2
()1
1
Fjω
ω
=
+
的傅里叶反变换f(t)
MATLAB程序如下:
symstw%定义两个符号变量t,w
Fw=sym('1/(1+w^2)');%定义频谱函数F(jw)
ft=ifourier(Fw,w,t);%对频谱函数F(jw)进行傅氏反变换
写出运行结果。
运行结果参考:
ft=
1/2*exp(-t)*Heaviside(t)+1/2*exp(t)*Heaviside(-t)
2、傅里叶变换的数值计算实现法
严格说来,如果不使用symbolic工具箱,是不能分析连续时间信号的。
采用数值计算方法实现连续时
间信号的傅里叶变换,实质上只是借助于MATLAB的强大数值计算功能,特别是其强大的矩阵运算能力
而进行的一种近似计算。
傅里叶变换的数值计算实现法的原理如下:
信号处理基础实验指导书——实验部分
8
对于连续时间信号f(t),其傅里叶变换为:
F(jω)
0
()jtlim()jn
n
fteωdtfneωτ
τ
ττ
∞∞−−
−∞→
=−∞
=∫=Σ
其中τ为取样间隔,如果f(t)是时限信号,或者当|t|大于某个给定值时,f(t)的值已经衰减得很厉害,可
以近似地看成是时限信号,则上式中的n取值就是有限的,假定为N,有:
F(jω)
1
0
()
N
jn
n
τfnτeωτ
−
−
=
=Σ
若对频率变量ω进行取样,得:
()()kFk=Fjω
1
0
()k0
N
jn
n
τfnτeωτkM
−
−
=
=Σ<<
通常取:
02
kkk
MM
ωπ
ω
τ
==,其中0
ω是要取的频率范围,或信号的频带宽度。
采用MATLAB实现上式
时,其要点是要生成f(t)的N个样本值f(nτ)的向量,以及向量e−jωknτ,两向量的内积(即两矩阵的乘积),
结果即完成上式的傅里叶变换的数值计算。
注意:
时间取样间隔τ的确定,其依据是τ必须小于奈奎斯特(Nyquist)取样间隔。
如果f(t)不是严格的
带限信号,则可以根据实际计算的精度要求来确定一个适当的频率0
ω为信号的带宽。
2-3用数值计算法实现上面门函数f(t)=ε(t+1)−ε(t−1)的傅里叶变换,并画出幅度频谱图.
分析:
该信号的频谱为F(jω)=2Sa(ω),其第一个过零点频率为π,一般将此频率认为是信号的带
宽。
但考虑到F(jω)的形状(为抽样函数),假如将精度提高到该值的50倍,即取05050Bω=ω=π,
则据此确定的Nyquist取样间隔为:
00
110.02
22
2
F
τ
ω
π
≤==
×
。
MATLAB程序如下:
R=0.02;%取样间隔τ=0.02
t=-2:
R:
2;%t为从-2到2,间隔为0.02的行向量,有201个样本点
ft=[zeros(1,50),ones(1,101),zeros(1,50)];%产生f(t