【精品】系统的复频域分析

1、实验六系统的复频域分析§ 6. 1 matlab函数lsim (用于系统函数)目的用lsim仿真由系统函数表征的因果lti系统的输出。基本题1. 定义系数向量al和bl用以描述由下面系统函数表征的因果lti系统：y-22. 定义系数向量a2和b2用以描述由下而系统函数表征的因果lti系统:3弘(s)二一5 + 0.33. 定义系数向量a3和b3用以描述由下面系统函数表征的因果lti系统：5 + 0.84. 利lsim和前面部分定义的向量求这些因果lti系统对由t=0:0.1:0.5,x=cos(t) 给出的输入的输出。以上四题解:匸0:0.1:0.5;x=cos(t);bl=l -2

2、;al=12;b2=3;a2=l 0.3b3=2a3=l 0.8yl=lsim(bl,al,x,t);subplot(2,2,l)plot(t,yl);y2=lsim(b2,a2,x,t);plot(t,y2); y3=lsim(b3,a3,x,t); subplot(2,2,3) plot(t,y3);§6.2作连续时间的零极点图目的这一节要学习如何在一个零极点图上展现有理系统函数的零极点。基本题1. 下列每个系统函数都对应于稳定的lti系统。用roots求每个系统函数的零极 点，如上所示的利用plot画出零极点图并作适当标注。(i) h(g 2曲s + 2s + 3z- *rr

3、/ 22+53 + 12(11) h(s) = s + 2$ + 10(iii)h(s) =2f +5s + 12u 2 +2s + l()($ + 2)若干不同的信号能有相同的拉普拉斯变换有理表达式，但有不同的收敛域。例如，具有单位冲激响应为hc(z) = eatu(t),hac(/) = -eatu(rt)的因果和反因果lti系 统就有相同的分子和分母多项式的有理系统函数为hc (5) = re(5)> -as + ahac =re(5)< -as + a然而，它们有不同的系统函数，因为它们有不同的收敛域。解：代码如下bl=15j;al=l 2 3;zsl=roots(bl)p

4、sl=roots(al)subplot(3,2,l)plot(real(zs 1 ),imag(zs 1), 'o*);gridaxis(卜 10 10-10 10)subplot(3,2,2)plot(real(psl ),imag(pslgridaxis(-10 10-10 10)b2=2 5 12;al=l 2 10;zs2=roots(b2)ps2=roots(a2)subplot(3,2,3)plot(real(zs2),imag(zs2),'o,);gridaxis(-10 10-10 10)plot(real(ps2),imag(ps2),'x')

5、;gridaxis(-10 10-10 10) b3=l 5; a3=l 2 3;zs3=roots(b3) ps3=roots(a3)subplot(3,2,5) plot(real(zs3)jmag(zs3)/o,); gridaxis(-10 10-10 10) subplot(3,2,6) plot(real(ps3)jmag(ps3)/xf);grid axis(-10 10-10 10)figure 1|复|目10file edit view insert tools desktop window help ns-1o1o q0-50510x100-5051010d-50510-1

6、0w2. 对1屮每个有理表达式，确定它们的收敛域。解：代码如下i.re(s)>-l;ii. re(s)为任意值 i.re(s)>-23. 对输入和输出满足下面微分方程：皱_3曲)=獰+ 2如+ 5口) dtdt2dt的因果lti系统，求系统函数的零点和极点，并完成一幅适当标注的零极点图。解：代码如下b=l -3;a=l 2 5;zs=roots(b)ps=roots(a)subpl ot(2,1,1)plot(real(zs),imag(zs),'o,);gridaxis(-10 10-10 10)subplot(2,l,2)plot(real(ps),imag(ps),&

7、#39;x');gridaxis(-10 10-10 10)中等题下面的练习需要用到函数plotpz,其具体实现如下。函数plotpz是用于对分 子和分母多项式的系数存入向量a和b中的lti系统画岀零极点图用的。这个函 数除了画图以外，还产生零极点的值。function ps,zs=plotpz(b,a，roc)ps=roots(a);zs=roots(b);ps=ps(:);zs=zs(:);maxl=max(abs(imag(ps;zs;j);maxr=max(abs(imag(ps;zs;i);plot(1.5*f-maxr maxrj0 0,w);hold ontext( 1.

8、5 * maxr,0/re');plot(0 0, 1.5*-maxi maxl,w);text(0,1.5 *maxl,tm)；plot(real(zs),imag(zs)/bo,);plot(real(ps) jmag(ps)/bx<);if (nargin('plotpz,)>2)% nargin是用来判断输入变量个数的函数if(any(real(ps)<roc)%检测矩阵中是否有非零元素,如果有，则返回1, 否则，返冋0ipole=max(real(ps(real(ps)<roc);plot(ipole ipole,1.5*-maxl maxlj

9、/r');endif(any(real(ps)>roc)rpole=min(real(ps(real(ps)>roc);plot(rpole rpole,1.5*-maxl maxl/r');endtext(roc,-l .25*maxl；roc');endaxis(,equar);grid;hold on任选宗量roc用来在图屮指岀收敛域。在该系统的收敛域内用选取roc是一个 点，plotpz将适当地给系统的收敛域标注。例如，试试执行下列命令»b=l -1;» a=l 3 1;» fps,zs=plotpz(b,a, 1);&

10、#187; ps,zs=plotpz(b,a,-2);b=l -u；a=l 3 1；ps,zs=plotpz(b,a,l);ps,zs=plotpz(b,a,-2);figure 1-n xfile edit viewinsert tools desktop windowhelp*u 口 3。1 n 1 代 a ® mg y 1因丄i ra4说明plotpz如何在roc内l田i出一个单一的点确定一个有理变换的收敛域。 解：通过any(real(ps)<roc)函数来确定极点与roc的关系，通过 min(real(ps(real(ps)>roc)函数找出边界，然后用plot

11、画出收敛域§6.3 matlab 函数 freqz目的学习用函数freqz绘制lti系统的频率响应。基本题1. 利用fieqz画出微分方程3心二+ 4如+曲)=尤字+ 5x描述的一个因果、 dt2dtdt2稳定的lti系统的频率响应的幅值和相位。解：代码如下：a=3 4 1;b=l 0 5;freqz(b,a);402000.40.50.60.80.9normalized frequency (xti rad/sample)mp)cdpwgew2%(saamap) as£dfigure 1亠口x! file edit viewinsert tools desktop win

12、dow helpd已m fe %兔的遐/复目旦60§6.4系统的时域和频域特性目的考虑几个出线性常系数微分方程描述的稳定的lti系统，对这些系统要求计 算它们的单位冲激响应和频率响应。尽管用频率响应或单位脉冲响应都足以完全 表征一个lti系统，但是将会明白有时候既从时域，又从频域来考虑系统特性是 很有利的。基本题 考虑由下面微分方程给出的一类因果lti系统dy(f)其屮d。 0以保证稳定性。定义系统i是满足上式，= 3的系统，系统ii是a0 = 1/3 的系统。1. 用解析法导出对应于上式的稳定lti系统的频率响应，同时确定这个频率响 应的幅值和相位。解：代码如下：al=l 3;bl

13、=3; freqs(bl,al) a2=l 1/3;b2=l/3;freqs(b2,a2)10frequency (rad/s)10(saa6ap) as?d0108o o o-24601010apn-e6e 乏-2-210ofrequency (rad/s)10'1 10°frequency (rad/s)101(saa6ap) ase-xzd0-2 o 1 o1010°10* =3 二旦0 51.522.53normalized frequency (xn rad/sample)4 2 0-2 eapnuew2. 定义w=linspace(0,10),利用fre

14、qz计算系统i和系统ii在w频率上的频率响 应，在单一的图上画出这两个频率响应的幅值。这两个幅值图与解析表达式屮频 率响应的幅值一致吗？解：代码如下：w=linspace(0,10);a=l 3b=3;freqz(b,a,w)丿 figure 1 xfile edit view insert tools desktop window helpq jj丄 a % % t© 3 x 宜(s8e<) ajd0.511.522.5normalized frequency (xn rad/sample)w=linspace(0,10); a=1 1/3; b=1/3; freqz(b,a

15、,w)6 8 0 2- - 1 1 - mp)詈uewo(sdajfap) as 启 do1020 0.511.522 53normalized frequency (xn rad/sample)20100结果分析：这两个幅值图与解析表达式中频率响应的幅值一致3用函数impulse计算系统i和系统ii在向量t=linspace(0,5)所定义的时间样本 点上的单位冲激响应。解：t=linspace(0,5);b=3;a=l 3j; impulse(b,a,t);t=linspace(0,5);b=l/3;a=l 1/3;impulse(b,a,t);file edit viewsmbinser

16、t tools desktop window help3a曰日©毀 ®跟謠=)旦4 单位冲激响应随时间衰减的速率与频率响应幅值随频率下降的速率之间是什 么关系？ ctft屮的何种性质说明这一关系？解：单位冲激响应随时间衰减的速率与频率响应幅值随频率下降的速率z间是一 对傅立叶变换关系中等题5. 巴特沃兹滤波器就属丁这样一类连续时间频率选择滤波器，它们能用线性常 系数微分方程实现。为了确定一个二阶畋=3的巴特沃兹滤波器的系数，键入 » wc=3;» b2,a2=butter(2,wc/s');宗量空表明所要求的是一个连续时间滤波器而不是离散时间滤波

17、器；如有缺损， 则由butter产生的是一个离散时间滤波器。向量b2和a2包含的是二阶微分方程 的系数，其格式与2.3节所用相同。用freqz计算并画出相应lti系统在 w=linspace(0,10)频率点上频率响应的幅值。在一张图上要包括2中已计算出的系 统i的频率响应的幅值。解：代码如下：wc=3;b2,a2=butter(2,wc,1s1); w=linspace (0,10);h2 = freqs(b2z a2,w)al=l 3;bl = 3;hl=freqs(bl,al,w)plot(w,abs (h2)r 1r *);hold on6. 利用频率响应幅值，哪个系统更逼近于哄=3的

18、理想低通滤波器，是系统i还 是二阶巴特沃兹滤波器？将每一系统的频率响应相位与该理想低通滤波器的频 率响应和位比较，情况怎样？ 解：二阶巴特沃兹滤波器更逼近于畋=3的理想低通滤波器；将每一系统的频率 响应相位与该理想低通滤波器的频率响应和位比较，相位基本保持-个固定值不 变深入题7. 用impulse计算二阶巴特沃兹滤波器在时间样本t=linspace(0,5)±的单位冲激 响应，在一张图上画出系统i和二阶巴特沃兹滤波器的单位冲激响应。解：代码如下：al=l 3;bl=3;t=linspace(0,5)h 1 =impulse(b l,al,t);wc=3;b2,a2=butter(2

19、, wc, 's*); h2=impulse(b2,a2,t);plot(t,hl,t')hold on;plot(t,h2);§6.5二阶系统的极点位置目的考查二阶系统阻尼比和无阻尼自然频率对系统零极点分布以及频率响应的 影响。基本题考查纠保持为1,:取不同值时的频率响应和极点位置。1.定义7g)百分别是而 = 0,1/44和2时,由上式确定的系统函数。 定义向量ala4分别是h|g)£($)分母多项式的系数向量。求出并画出每个 系统极点的位置。解：代码如下：al=l 0 1;a2=l 1/2 1;a3=l 2 11;a4=l 4 1j;psl=roots

20、(al);ps2=roots(a2);ps3=roots(a3);ps4=roots(a4);subplot(4,l,l);plot(real(al ),imag(a 1 ),x);subplot(4,l,2);plot(real(a2),imag(a2),'x,);subpl ot(4,1,3)；pl ot(real(a3),i mag(a3 ),'x');subpl ot(4,1,4);plot(real(a4),imag(a4)/xr);llllilllt111it c00_10.20.30.405060.70.80.91-111ii t1ii11i111111t

21、.50 550 60 650.70.750 80 850.90 9511111t11.21.31.41.51.61.71.81.92iiiiii0 : 1t :彳 iiiiii11.522.533.542. 泄义omega二卜5:0.1:5是这些频率，在这些频率上要计算这4个系统的频率响 应用freqz函数计算并画出在1屮所定义的这4个系统的|丹伽)|。滞 g和01 的频率响应定性来看是怎样的不同？你能解释系统的极点位置为何会导致这种 差别？另外，从儿何上你能论证为什么对所有这4种系统都有相同的片(.3)爲值 吗？解：代码如下：al=l 0 1;a2=ll 1/2 1;a3=l 2 1;a4=

22、l 4 1;b=l;omega=-5:0.1:5;h 1 =abs(freqs(b,a 1,omega);h2=abs(freqs(b,a2,omega);h3=abs(freqs(b,a3,omega);h4=abs(freqs(b,a4,omega);subplot(4,l,l);plot(hr);subplot(4,l,2);plot(h2);subplot(4,l,3)；plot(h3);subplot(4,l,4);plot(h4);结果分析:中等题八2下而将耍画出hac(s) = -当改变：和时极点位置的轨迹，并 s + s + ct)n且看看这些参数变化是如何影响系统的频率响应的

23、。3. 首先保持® =1,在()5：m1()范围内改变：。定 zetarange=o logspace( 1,1,99) 以得到在0<c <10内100个按对数分开的点。定义azata是一个3x10 0的矩阵， 其中每一列是当：为zetarange中对应列的值时，h(s)分母多项式的系数。定义 zetapoles是2x100的矩阵，其中每一列就是azeta对应列的根。在单独的图上， ifflitb zetspoles每一行实部对于虚部的图，并叙述i田i图的轨迹。在图上指出下列 各点：$ = 0,1/4和2。为了在图中得到具有等长度轴的平方纵横比，可以键入 » a

24、xis(,equal,);» axis(-4 0 -2 2);定性叙述一下，当：从0改变到1,然后又从1改变到10时，估计频率响应会如 何变化。解：代码如下：wnhrrnulop zaarangeuoogspace(l1999)k azaehones(3 -100)八 azsap-hze erange - zsapo 一 eshzerosp100)八 for nurrn 八 zeapooscd)hroois(azatacb)- end八 subplot(211) p 一2(real(zesp2es( 1o) jmag(zaap2es( 1o)；axis(t0 2 2)八 subplo

25、t(2?2)八 ploi(r<l(zetapoles(2o)jmag(zetapoos(2o)；axis(l-4 0 ,2 2r4 .保持1/4而要考查®，从0增加到10的效果。定义omegarange二0 logspace(-l,l,99)是按对数分开的在所关心的区域的100个点。用类似于在3中 定义a zeta和zetapoles的方式定义aomega和omegapoles。在一张单独的图上， 对omegapoles的每一列画出实部对虚部的图。你如何估计变化®是怎样改变频 率响应h(je)的？利用freqz求当® =2和"1/4时的频率响应，并

26、画出这个频率 响应的幅值。将这个结果与2屮对纬=1和"1/4的图进行比较。它们有什么不 同？这个与从按omegepoles画出的轨迹图所估计的符合吗？解：代码如下：seta=l/4;n=100;omegarange=fo logspace(-1,1,99);aomega=ones(3j 00);aomega(2,: )=1/2. *o megarange;aomega(3, :)=omegarange. a2;omegapoles=zeros(2,100);for n=l:n;omegapoles(: ,n)=root s(aomega(:, n);end;subplot(4,l,l

27、);plot(real(omegapoles(l/)5imag(omegapoles(l,:);subplot(4,l,2);plot(real(omegapoles(2,:)jmag(omegapoles(2,:);wn=2;seta=l/4;b=4;a=l 1 4;fi*eqs(b,a,omegarange);10'10°10'110-210°frequency (rad/s)100101(s<da>6a>p)a>s?do o 1 20o o o5 0 5- 1 1深入题5. 在上式屮的：没有理由一定要是止的。对于位于1o和0之间

28、的：重做3。当：是负的时，由/描述的系统能够既是因果的又是稳定的吗？另外，用freqz画!ll = -l/4和时该系统频率响应的幅值。具有由freqz计算出的频率响应的 这个系统是因果的吗？另外，解释一下这个图与在2中对 = 1/4画出的频率响应 的幅值之间的任何异同点。解：代码如下：wn=l;n=100;zetarange=-o logspace(-l ,1,99);azata=ones(3 j 00);azata( 2, :)=zetarange;zetapoles=zeros(2,100);for n=l:n;zetapoles(:,n)=roots(azata(:,n);end;sub

29、plot(2,l,l);plot(real(zetapoles( 1,:) ,imag(zetapoles(l,:);subplot(2,l,2);plot(real(zetapoles(2,:)5imag(zetapoles(2,:);seta=-l/4;a=l -1/2 1j;b=l;freqs(b,a,zetarange);402%-20246810normalized frequency (xz rad/sample)20oo20(sq>(d6<dp)(dsplld00-3-2o ooo o325-2-1.5-1050normalized frequency (xz rad

30、/sample)§ 6. 6实现非因果连续时间滤波器目的讨论由因果和反因果滤波器实现非因果滤波器。i.如果心(s)是与(6.4)式有关的系统函数，比(s)是与(6.3)式有关的系统函数,这两个系统的极点存在什么关系？比和hac(s)是什么关系？解：两个系统的极点一个取最大值，一个取最小值，并且方向相反，hc(s)=-hac(s)； 对于下面的题，将考虑微分方程+2y(/) = x(/)(6.5)dt2.求与(65)式有关的系统函数h(s)和全部可能的h(s)的收敛域，对每一种收敛域，确定相应的lti系统的单位冲激响应。解：代码如下a=l 2;b=l;t=linspace(0,10);

31、yl 二 impulse(b，a,t);y2=impulse(-b,a,t);subplot(2,l,l);plot(yl);subplot(2,l,2);plot(y2);当系统为因果系统时，h(s)的收敛域为(2, 00),冲击响应为e2lu(t)；当系统为反因果系统时,h(s)的收敛域为(8,2),冲击响应为e2tu(-t)；3. 对2中的每一种收敛域，对应的微分方程的辅助条件是什么？解:辅助条件分別是因果系统、非因果系统4定乂 a和b包含了 h($)分子和分母多项式的系数。对于具有系统函数为h(s)的 因果系统，用impulse证实单位冲激响应的解析表达式。将有impulse产生的在

32、时间样本t=-5:0.01:5上的单位冲激响应的值存入向量h中。画出h对t的图。应注意，因为hs=impulse(b,a,ts)产生的是对罚-拾)响应的样木，所以将不得不 对impulse的输入选择适当的时间样本，并且更将结果补上适当个数的零。解：代码如下：b=l;a=l 2;t=0:0.01:5;%k=length(t);%ts=t zeros(l,k);y=impulse(b,a,t)；plot(y);5.对反因果系统重做4。记住，impulse(b,a,ts)是假设在a和b中的系数对应于 一个因果系统。需要对一个时间倒置系统定义一组新的系数，然后将由impulse 计算出的单位冲激响应适

33、当地翻转。解：代码如b=l;a=-l 2j;ts=o:o.ol :5;t=-5:0.01:5;h=zeros(l,1001);p=impulse(b,a,ts);hs=conj(p');for n=l :501;h(n)=hs(502-n);end;plot(t,h);6. 用解析方法计算当输入曲)"“2心)时，满足(6.4)式的反因lti系统的输出。解：代码如下：b=l;a=-l 2;ts=0:0.01:5;t=-5:0.01:0;h=zeros(l,1001);p=impulse(b,a,ts);hs=conj(p*);for n=l :501;h(n)=hs(502-n

34、);end;x=exp(5 *t./2).*heaviside(-t);y=conv(x,h);7. 用lsim验证在6中导得的反因果系统在时间样本t上的输出。和impulse 样，函数lsim(b,a,x,ts)也是假设向量n和b对应于一个因果系统，所以必须用时 间倒置的微分方程来做。这个微分方程的系数应该在5屮己计算出。注意，时间 倒置系统的输入也必须是时间倒置的。解：代码如下：b=l;a=-l 2j;ts=-5:0.01:0;=5:0.01:5;x=exp(5 *ts./2). *heaviside(-ts);h=zeros(l,1001);hs=lsim(b,a,x,ts);for n

35、=l:501;h(n)=hs(502-n);end;plot(t,h);中等题现考虑三阶微分方程獰 + 24 2-26,(0 =牌 + 7 如 + 21x(/)(6.6)dt3 dt2dt dt2dt8. 求与(6.6)式有关的系统两数h(s),并按6.2节画出这个系统两数的零极点图。确定全部可能的收敛域。对于哪一种收敛域该系统是稳定的？解：代码如下：a=l 1 24 -26;b=l 721;zs=roots(b);ps=roots(a);t=linspace(0j0);hl=impulse(b,a,t)plot(real(zs),imag(zs),o,);hold onplot(real(p

36、s),imag(ps),'x');grid结果分析：对于resvl时，该系统是稳定的9. 用residue确定h(s)的部分分式展开式。对由8所确定的每个收敛域，用解析法确泄有关的单位冲激响应。解：代码如卜：a=l 1 24 -26;b=l 721;h,p,k=residue(b,a)10. 对确(6.6)式的因果系统，用impulse验证曲9计算出的相应的单位冲激响 应，时间样本用向量t。与这个因果系统有关的曲)上的辅助条件是什么？ 解：代码如下：a=l 1 24 -26;b=l 721;t=0:0.01:10;hl=impulse(b,a,t)plot(hl);11对满足(

37、6.6)式的反因果系统，用impulse验证由9计算出的相应的单位冲激 响应，时间样本用向量t。记住，必须首先确立这个时间倒置的微分方程的系数,然后从这个时间倒置方程的因杲单位冲激响应导出(6.6)式这个反因杲的单位冲激响应。与这个反因果系统有关的y上的辅助条件是什么？解：代码如下：a=-l 1-24 -26;b=l -7 21;ts=0:0.01:5;t=0:0.01:10;h=zeros(l,1001);p=impulse(b,a,ts);hs=conj(p');for n= 1:501;h(n)=hs(502-n);end;plot(t, h);grid;深入题12 并行实现一个lti系统是很容易这样实现的：将系统函数分解为h(s) = h(s) + h2(s),以使得系统输出w)能分别由计算对hl和丹2($)的响应而计算出來。注意，7和丹2($