MATLAB讲义第二章.docx
- 文档编号:9382252
- 上传时间:2023-02-04
- 格式:DOCX
- 页数:19
- 大小:201.30KB
MATLAB讲义第二章.docx
《MATLAB讲义第二章.docx》由会员分享,可在线阅读,更多相关《MATLAB讲义第二章.docx(19页珍藏版)》请在冰豆网上搜索。
MATLAB讲义第二章
第二章 序列的付氏变换(DTFT)和Z变换
一、序列的DTFT:
X(ejw)
序列x(n)的付氏变换如下:
方法一:
定义法此法只为让大家熟悉DTFT形成的过程
因为MATLAB无法计算连续变量w,只能在
(或考虑对称性,在
)范围内,把w赋值为很密的、长度很长的向量来近似连续变量。
通常最简单的就是赋以K个等间隔的值:
它表示把基本数字频率的范围
均分成K份后,每一份的大小。
则频率向量表示为w=[w1,w2,…,wk]=k
dw,
序列的位置向量为n=[n1,n2,…,nN]
则DTFT的计算式可用一个向量与矩阵相乘的运算来实现:
[X(w1),X(w2),…,X(wk)]
=[x(n1),x(n2),…,x(nN)]
=[x(n1),x(n2),…,x(nN)]
又因为(jnT*w)=jnT*kdw=j*dw*n’*k
所以
Xw=xn*exp(j*dw*n’*k)
例1:
求有限长序列xn=[13531],n=-1:
3的DTFT,画出它在w=-8~8rad/s范围内的频率特性,讨论其对称性。
再把xn左右移动,讨论时移对DTFT的影响。
解:
根据DTFT定义得
将w在-8~8rad/s之间分为1000份。
>>xn=[13531];nx=-1:
3;
>>w=linspace(-8,8,1000);%设定频率向量¿
>>Xw=xn*exp(-j*nx'*w);%用定义计算DTFT
>>subplot(3,5,1);stem(nx,xn);axis([-2,6,0,6]);title('ÔʼÐòÁÐ');ylabel('xn');%画序列图
>>subplot(3,5,2);plot(w,abs(Xw));title('·ù¶È');%画幅频曲线
>>subplot(3,5,3);plot(w,angle(Xw));title('Ïà½Ç')%画相频曲线
>>subplot(3,5,4);plot(w,real(Xw));title('ʵ²¿')%画实频曲线
>>subplot(3,5,5);plot(w,imag(Xw));title('Ð鲿')%画虚频曲线
>>nx_2=nx+2;%使xn右移2位
>>Xw_2=xn*exp(-j*nx_2'*w);
>>subplot(3,5,6);stem(nx_2,xn);axis([-2,6,0,6]);ylabel('xn_2');title('ÓÒÒÆ2λ')
>>subplot(3,5,7);plot(w,abs(Xw_2));
>>subplot(3,5,8);plot(w,angle(Xw_2));
>>subplot(3,5,9);plot(w,real(Xw_2));
>>subplot(3,5,10);plot(w,imag(Xw_2));
>>nx_1=nx-1;%使xn左移1位
>>Xw_1=xn*exp(-j*nx_1'*w);
>>subplot(3,5,11);stem(nx_1,xn);axis([-2,6,0,6]);ylabel('xn_1');title('×óÒÆ1λ')
>>subplot(3,5,12);plot(w,abs(Xw_1));
>>subplot(3,5,13);plot(w,angle(Xw_1));
>>subplot(3,5,14);plot(w,real(Xw_1));
>>subplot(3,5,15);plot(w,imag(Xw_1));
由上述结果得出以下结论:
(1)序列的DTFT是连续函数。
(2)序列的DTFT是周期函数,周期为2pi。
(3)本题给出的是实序列,实序列的DTFT具有对称性,幅频和实频偶对称;相频和虚频奇对称。
(4)信号的时移不影响幅频特性,只影响相频。
方法二:
调用内部函数法(Freqz)
MATLAB工具箱中有计算DTFT的专用函数freqz。
它的调用方式为:
(1)[h,w]=freqz(b,a,N)
h—复频率响应
N—N点复频率响应,是把
分割的份数,计算的是频率部分的特性。
若省略N,默认512
b,a—滤波器系数的分子分母向量,且以z-1的升幂排列的系数向量b为输入序列xn的系数向量,a为输出序列yn的系数向量,且a0要归一化,即a0=1.
w—数字频率向量,它把0-π均分为N份,w=[0:
N-1]π/N
此式可计算出正频率区间特性。
(2)[h,w]=freqz(b,a,N,’whole’)
用于画出全奈氏频率范围内的特性。
此时N把0-2π均分。
(3)[h]=freqz(b,a,w)
在自己选定的频点向量W上计算频率特性。
(4)freqz(b,a)
若没有给出左端,则直接在当前窗口给出频率响应的幅频响应和相频响应。
例:
差分方程:
y(n)-0.9y(n-1)=0.5x(n)+0.8x(n-1),求它的频率响应H(ejw),并画图。
并求输入为x(n)=cos(0.1πn)u(n)的稳态输出ys.
>>b=[0.5,0.8];
>>a=[1,-0.9];
>>[h,w]=freqz(b,a);%求出频率响应(0到π分成500点)
>>subplot(2,1,1),plot(w,abs(h)),gridon;title('|h(ejw)|')
>>subplot(2,1,2),plot(w,angle(h)),gridon;title('angle[h(ejw)]')
>>
稳态输入x(n)的频率为wx=0.1π,初始相角为θx=0,系统在该频点处的响应为:
>>wx=0.1*pi;
>>nwx=floor(wx/pi*512)
nwx=
51
>>hwx=h(nwx)
hwx=
1.2085-4.0139i
>>Ax=abs(hwx)
Ax=
4.1919
>>thetax=angle(hwx)
thetax=
-1.2784
>>y=Ax*cos(0.1*pi*n-thetax);
>>subplot(2,1,1),stem(n,x),title('x');
>>subplot(2,1,2),stem(n,y),title('y');
>>subplot(2,1,1),stem(n,x,'*'),title('x');
>>subplot(2,1,2),stem(n,y,'.'),title('y');
由运行结果Ax=4.1919,thetax=-1.2784得
即
则ys=4.2694cos(0.1πn-1.2784)=4.1919cos[0.1π(n-4.0693)]
三、Z变换的计算
(1)两个序列的Z变换
例1:
x1(n)=[2,3,4],ns1=0;x2(n)=[3456],ns2=0;求两序列的Z变换,及两序列的卷积y(n)及Y(z).
由Z变换的定义
知,
X1(z)=2z0+3z-1+4z-2
X2(z)=3z0+4z-1+5z-2+6z-3
∵y(n)=x1(n)*x2(n)
∴Y(z)=X1(z)X2(z)
∴Y(z)=(2z0+3z-1+4z-2)(3z0+4z-1+5z-2+6z-3)
=6+17z-1+34z-2+43z-3+38z-4+24z-5
则,y(n)=[61734433824]
y(n)=x1(n)*x2(n)=conv(x1,x2)
求卷积程序:
>>x1=[234];nx1=0:
2;x2=[3456];nx2=0:
3;
>>y=conv(x1,x2)
y=
61734433824
即y(n)=[61734433824]
从而Y(z)==6+17z-1+34z-2+43z-3+38z-4+24z-5
由以上看出,两序列的卷积和它们的Z变换的相乘结果是一致的。
(时域卷积定理)
因为Z变换的计算还涉及到收敛域的确定,这是较为复杂的,所以目前没有一个计算Z变换的很好的函数.
三、Z反变换
MATLAB的内部函数residuez可计算离散系统的极点留数。
调用格式:
[r,p,C]=residuez(b,a).
b,a――H(z)中,分子分母多项式以z-1的升幂排列的系数向量。
p――极点向量,即分母的根
r――部分分式分解后,各极点所对应的留数Ak
C――无穷项多项式系数向量,仅在M≥N时存在。
则
则Z反变换为:
y(n)=r
(1)p
(1)nu(n)+…+r(N)p(N)nu(n)+C
(1)δ(n)+C
(2)δ(n-1)+…
若H(z)是系统函数,则其反变换h(n)就是它的单位脉冲响应,可用函数impz来计算,可用该命令来核对计算。
其调用方式为:
h=impz(b,a,L)
其中L——待求脉冲序列h(n)的长度
实例1:
计算的
反变换
先用函数poly求出分母多项式的系数,再用residuez求X(z)的极点和留数
程序:
>>b=1;a=poly([0.9,0.9,-0.7]);
>>[r,p,C]=residuez(b,a)
则其z变换为:
则
r=
0.2461
0.5625
0.1914
p=
0.9000
0.9000
-0.7000
C=
[]
四、用Z变换求系统输出
(1)方法1:
采用residuez函数
实例1:
系统函数
,输入信号x=[2345],用Z变换求出输出y(n)。
解:
∵Y(z)=H(z)·X(z)=
=
∴ B(z)是X(z)与-3z-1的卷积,用函数conv可求得
求出Y(z)后,用函数residuez可求其反变换,即得y(n)。
>>x=[2345];nsx=0;
>>nfx=length(x)-1;
>>b=-3;a=[2-52];
>>B=conv(-3,x);A=a;
>>[r,p,k]=residuez(B,A)
r=
-10.2500
32.0000
求得
则
p=
2.0000
0.5000
k=
-24.7500-7.5000
>>nf=input('nf');
nf8
求得Yi(z)=
Yd(z)=-24.75-7.5z-1
>>n=0:
nf;
>>yi=(r
(1)*p
(1).^n+r
(2)*p
(2).^n).*stepseq(0,0,nf);
>>yd=k
(1)*impseq(0,0,nf)+k
(2)*impseq(1,0,nf);
>>y=yi+yd
>>subplot(1,2,1),stem([nsx:
nfx],x);title('x'),line([0,10],[0,0]);
>>subplot(1,2,2),stem([0:
nf],y);title('y'),line([0,10],[0,0]);
(2)方法2:
采用filter函数
若不需求出表达式,可用函数filter来解。
y=filter(b,a,x)
其中:
b,a——差分方程中x(n),y(n)的系数数组
b=[b0,b1,…,bM]
a=[a0,a1,…,aN]
x——输入序列
当初始条件为0,输入序列为
时,y即为系统的单位脉冲响应h(n)
例:
>>x=[2,3,4,5,zeros(1,10)];b=-3;a=[2-52];
>>y=filter(b,a,x);
>>ny=length(y)-1
ny=
13
>>subplot(1,2,1),stem([0:
ny],x);title('x'),line([0,10],[0,0]);
>>subplot(1,2,2),stem([0:
ny],y);title('y'),line([0,10],[0,0]);
>>
五、Z变换解差分方程
1、递推法(单位脉冲响应)
法1:
采用filter函数
已知a,b向量,给定输入x(n),则可求得输出y(n),其长度与输入x(n)相同。
当初始条件为0,输入为单位脉冲函数时,系统的输出为脉冲响应h(n)。
>>b=2*[1,2,3];
>>a=[1,-0.8,0.4];
>>Y=[-2,-3];
>>X=[1,1]
X=
11
>>clear
>>b=1;
>>a=[1,-1,0.9];
>>x=[1,zeros(1,7)];
>>y=filter(b,a,x)
y=
1.00001.00000.1000-0.8000-0.8900-0.17000.63100.7840
法2:
采用impz函数
调用方法:
(h,t)=impz(b,a):
返回参数h是脉冲响应的数据,t是脉冲响应的采样时间间隔
(h,t)=impz(b,a,N):
N为指定脉冲信号的长度。
(h,t)=impz(b,a,n,Fs):
Fs指定脉冲信号的采样频率,其默认值为1。
例:
>>x=[2,3,4,5,zeros(1,10)];b=-3;a=[2-52];
>>N=14;
>>[h,t]=impz(b,a,N);
h=
1.0e+004*
-0.0001
-0.0004
-0.0008
-0.0016
-0.0032
-0.0064
-0.0128
-0.0256
-0.0512
-0.1024
-0.2048
-0.4096
-0.8192
-1.6384
t=
0
1
2
3
4
5
6
7
8
9
10
11
12
13
>>stem(t,h),title('h(n)=impz(b,a,14)')
六、H(z)的零极点
roots函数求其零极点,
zplane(b,a)可绘出零极点图
例:
,求零极点,并作图。
>>b=[2040-2];
>>a=[10.3000];
>>z=roots(b)
极点
z=
-0.0000+1.5538i
-0.0000-1.5538i
0.6436
-0.6436
>>p=roots(a)
零点
p=
0
得到传输函数
0
0
-0.3000
>>zplane(b,a)
>>
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- MATLAB 讲义 第二