数字信号处理 第六章课后MATLAB答案 华工

第六章习题与思考题参考答案1. 简述飞机红外图像识别中用到的五个红外特征量各自的作用？1）长宽比：反应了目标的几何形状；2）复杂度：反应了红外目标轮廓的情况；3）紧凑度：反应了红外目标在其所在最小外接矩形中所占比重；4）均值对比度：反映了目标的物理特性与背景的物理特性之间的关系；5）部分最亮像素点数与目标总像素数的比值：反应了目标的明暗变化情况；2. 可视密码共享中,如果实现（4 , 5）门限的可视密码分享，程序将做如何编写？clear allclose allM=imread('0.jpg');ss=rgb2gray(M);figureimshow(ss);[m n]=size(ss);for i=1:m*nif ss(i)>250ss(i)=250;endends=double(ss)+1;x=[1 2 3 4 5];g1=zeros(m,n);g2=zeros(m,n);g3=zeros(m,n);g4=zeros(m,n);g5=zeros(m,n);yy1=zeros(m,n);yy2=zeros(m,n);yy3=zeros(m,n);yy4=zeros(m,n);yy5=zeros(m,n);y1=zeros(m,n);y2=zeros(m,n);y3=zeros(m,n);y4=zeros(m,n);y5=zeros(m,n);for j=1:m*na1=mod(2*j,251);a2=mod(3*j,251);a3=mod(5*j,251);f=[a1 a2 a3 s(j)];g1(j)=polyval(f,x(1));yy1(j)=mod(g1(j),251);g2(j)=polyval(f,x(2));yy2(j)=mod(g2(j),251);g3(j)=polyval(f,x(3));yy3(j)=mod(g3(j),251);g4(j)=polyval(f,x(4));yy4(j)=mod(g4(j),251);g5(j)=polyval(f,x(5));yy5(j)=mod(g5(j),251);endy1=uint8(yy1-1)y2=uint8(yy2-1);y3=uint8(yy3-1);y4=uint8(yy3-1);y5=uint8(yy3-1);figure,imshow(y1);figure,imshow(y2)figure,imshow(y3);figure,imshow(y4);figure,imshow(y5);l1=(x(2)*x(3)*x(4)*x(5))/[(x(1)-x(2))*(x(1)-x(3))*(x(1)-x(4))*(x(1)-x(5))];l2=(x(1)*x(3)*x(4)*x(5))/[(x(2)-x(1))*(x(2)-x(3))*(x(2)-x(4))*(x(2)-x(5))];l3=(x(1)*x(2)*x(4)*x(5))/[(x(3)-x(1))*(x(3)-x(2))*(x(3)-x(4))*(x(3)-x(5))];l4=(x(1)*x(2)*x(3)*x(5))/[(x(4)-x(1))*(x(4)-x(2))*(x(4)-x(3))*(x(4)-x(5))];l5=(x(1)*x(2)*x(3)*x(4))/[(x(5)-x(1))*(x(5)-x(2))*(x(5)-x(3))*(x(5)-x(4))];rr1=zeros(m,n);r=zeros(m,n);for j=1:m*nrr1(j)=mod(yy1(j)*l1+yy2(j)*l2+yy3(j)*l3+yy4(j)*l4+yy5(j)*l5,251);endr=uint8(rr1-1);figure,imshow(r);3. 已知图像⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡=805020016821024015010017018013070901262068M 。

《数字信号处理》第6章课后作业答案6.1（1）已知IIR 数字滤波器的系统函数为 (1) 232164016()81061z H z z z z -+=-+- 试写出滤波器的差分方程，并分别画出直接I 型、直接Ⅱ型、转置直接Ⅱ型、级联型和并联型结构图。

解:经化解原式可得：123123252()5311448z z z H z z z z -------+=-+-直接I 型：直接Ⅱ型:级联型：注意，对于级联型，一定要化成负幂次，再写系数！经对原式进行分解得：11211221 2.5()110.2512z z z H z z z z-------+=⨯--+并联型： 注意：系数b,a 是()H z z的系数！ b=[0,0,16,-40,16]; a=[8,-10,6,-1,0]; [K,z,d]=residue(b,a) KK1=[K(1),K(2)]; zz1=[z(1),z(2)];[b2,a2]=residue(KK1,zz1,0) 经原式分解得:111211.2 4.8 5.6()1610.2510.5z H z z z z -----=-++--+6.1（2）略 6.4Matlab 程序： clear; fp=5000; wp=2*pi*fp; fs=10000; ws=2*pi*fs; ap=3; as=30;[N,wc]=cheb1ord(wp,ws,ap,as,'s'); [B,A]=cheby1(N,ap,wc,'s') freqs(B,A);系统函数：17439213171.220510()18271 1.153910 1.25510 1.72410H s s s s s ⨯=++⨯+⨯+⨯ 图：6.6试设计一个巴特沃斯型模拟带通滤波器，并用Matlab 验证结果，要求带宽为200Hz ，中心频率为1000Hz ，通带内衰减不大于3dB,在频率小于830Hz 或大于1200Hz 处的衰减不小于25dB. 解：（1） 模拟带通滤波器的技术指标要求为：BW Ω=400πrad/s; 0Ω=2210πrad/s; p α=3dB; s α=25dB; 因为：ph pl Ω-Ω=400π；ph pl Ω•Ω=4000000π2； 所以可得：ph Ω=2210πrad/s; pl Ω=1810πrad/s; sh Ω=2400πrad/s; pl Ω=ll1810πrad/s; (2) 归一化频率为： 2210 5.525;400ph ph BWηΩ===Ω 18104.525;400pl pl BW ηΩ===Ω 24006;400sh sh BW ηΩ===Ω 1660 4.15;400sl sl BW ηΩ===Ω 2025;sh sl ηηη==（3）归一化原型模拟低通滤波器()a G p 的技术指标要求为：1;p λ= 222200min sh sl s sh sl ηηηηληη⎡⎤--=⎢⎥⎣⎦， ；s λ=1.83；（4）设计归一化原型低通滤波器：()a G p0.10.110.12101lg 101 4.762;2lg (101) 1.029;p s ss Ns N αααλλλ-⎛⎫- ⎪-⎝⎭≥-==-=c 所以 N=5;（5）查表得： 23451()1 3.3261 5.236 5.2361 3.2361an G u u u u u u=+++++ （6）()()ca an pu G p G u λ==（7）2()()()ph pl ph pla s p s Ha s G p Ω-Ω=+ΩΩ=Matlab 程序：fp=[905,1105]; fs=[830,1200]; wp=2*pi*fp; ws=2*pi*fs; ap=3; as=25;[N,wc]=buttord(wp,ws,ap,as,'s'); %巴特沃斯型模拟带通滤波器 [B,A]=butter(N,wc,'s'); f=500:1500; w=2*pi*f;H=freqs(B,A,w); subplot(2,1,1);plot(f,20*log10(abs(H))); grid on;axis([500,1500,-80,5]); xlabel('f/Hz'); ylabel('幅度/dB'); subplot(2,1,2); plot(f,angle(H)); grid on;axis([500,1500,-5,5]); xlabel('f/Hz'); ylabel('相位/dB');6.7 解：0.1T s =;112()(2)(3)23a s H s s s s s +-==+++++;122;3;s s ∴=-=-23122131(),,2()11T T T T H z z e z e T TH z e Z e Z ------==-∴=+--所以相应的的极点为Matlab 程序： clear; b=[1,1]; a=[1,5,6]; Fs=10;[B,A]=impinvar(b,a,Fs); [H,w]=freqz(B,A,'whole'); plot(w/pi,20*log10(abs(H)));6．8试用双线性变换法设计一个巴特沃斯型低通数字滤波器，并用matlab ，验证结果，给定技术指标为100,300,3,20,p s p s f Hz f Hz dB dB αα====采样频率为1000Hz 。

数字信号处理MATLAB习题数字信号处理MATLAB 习题M1-1 已知1()cos(6)g t t π=，2()cos(14)g t t π=，3()cos(26)g t t π=，以抽样频率10sam f Hz =对上述三个信号进行抽样。

在同一张图上画出1()g t ，2()g t 和3()g t 及抽样点，对所得结果进行讨论。

解：从以上两幅图中均可看出，三个余弦函数的周期虽然不同，但它们抽样后相应抽样点所对应的值都相同。

那么这样还原回原先的函数就变成相同的，实际上是不一样的。

这是抽样频率太小的原因，我们应该增大抽样频率才能真实还原。

如下图：f=50Hz程序代码f=10;t=-0.2:0.001:0.2;g1=cos(6.*pi.*t);g2=cos(14.*pi.*t);g3=cos(26.*pi.*t);k=-0.2:1/f:0.2;h1=cos(6.*pi.*k);h2=cos(14.*pi.*k);h3=cos(26.*pi.*k);% subplot(3,1,1);% plot(k,h1,'r.',t,g1,'r');% xlabel('t');% ylabel('g1(t)');% subplot(3,1,2);% plot(k,h2,'g.',t,g2,'g');% xlabel('t');% ylabel('g2(t)');% subplot(3,1,3);% plot(k,h3,'b.',t,g3,'b');% xlabel('t');% ylabel('g3(t)');plot(t,g1,'r',t,g2,'g',t,g3,'b',k,h1,'r.',k,h2,'g.',k,h3,'b.')xlabel('t');ylabel('g(t)');legend('g1(t)','g2(t)','g3(t)');M2-1 利用DFT的性质，编写一MATLAB程序，计算下列序列的循环卷积。

第六章 快速傅里叶变换(FFT)1. 如果一台通用计算机的速度为平均每次复乘需100μs,每次复加需20μs,今用来计算N=1024点的DFT[x(n)],问用直接运算需要多少时间,用FFT 运算需要多少时间。

解：ss FFT ss FFT N N a N N N m ss DFT s DFT DFT N a m FFT FFT DFT DFT μμμμμμμ23221027682666221020482010241051210512010240log ),10242(,5120log 210820104,10410104,104104102444⨯=⨯⨯=⨯======⨯=⨯⨯⨯=⨯⨯⨯⨯=⨯===作复加所需时间作复乘所需时间作复加所需时间作复乘所需时间作复乘 2. 用图6.8所示流程图验证图6.7所示的8点变址运算。

证明：由图6.8知取A=x(0),B=x(4)N=8X(k)=12/,,1,0),()(21-=+N k k X W k X k NX(N/2+k)=12/,,1,0),()(21-=-N k k X W k X k N5.试证实以下流图是一个N=8的FFT 流图.其输入是自然顺序的,而输出是码位倒置顺序的,试问这个流图是属与时间抽取法还是频率抽取法?并比较与书中哪一个流图等效。

解：这个流图属于频率抽取法。

6.试设计一个频率抽取的8点FFT 流图,需要输入是按码位倒置顺序而输出是按自然顺序的。

解:设计的流图为第五题的流图左右翻转180度。

∑∑-=-==+=1202/21202/1)()12()()2(N k kr N N k kr N W k x r X W k x r X7.试用图6.14(a)中的蝶形运算设计一个频率抽取的8点IFFT 流图。

解：X(0) 1/2 x(0) X(4) x(1)X(2) x(2)X(6) x(3)X(1) x(4)X(5) x(5)X(3) x(6)X(7) x(7)9.试作一个N=12点的FFT 流图,请按N=2,2,3分解,并问可能有几种形式?解：可能有三种先分成2组，每组有6各点，后每组内再分成两组322⨯⨯=N时间顺序为x(0),x(4),x(8),x(2),x(6),x(10),x(1),x(5),x(9),x(7),x(11)频域顺序为X(0),X(1),X(2),X(3),X(4),X(5),X(6),X(7),X(8),X(9),X(10),X(11)流图如图6.18解：由题可得∑∑-=-=-=-=∴-⋅⋅⋅====102210)(|)(1,,1,0,)()(N n kn Nj z z k N j k N n ne n x z X N k e z z z n x z X k ππ由于(a)将M 点序列分成若干段N 点序列，设段数为k 即N k M kN )1(-≥>并令kn N j N n k i i z z k en y z X N n N k M N k M n N k n x n y N n x n y n x n y k π21010110)]([[|)(11)1(,01)1(0],)1([)()()()()(--=-==-∑∑=⎩⎨⎧-≤<------≤≤-+=+==若用N 点FFT 计算)(k z X 先由x(n)形成)(n y i ，再计算∑-=10)(k i i n y 的N 点FFT 即可(b)先将序列添加一点等于零的点，使得⎩⎨⎧-≤≤-≤≤=1,010),()(0N n M M n n x n x再计算)(0n x 的N 点FFT 即10,)(|)(20-≤≤=∑-=N k e n x z X kn N j z z k π即可13.已知X(K),Y(K)是两个N 点实序列x(n),y(n)的DFT 值,今需要从X(K),Y(K)求x(n),y(n)值,为了提高运算效率试设计用一个N 点IFFT 运算一次完成。

第六章
6.2 如何进行一下操作：
（1）翻转模块
（2）给模型窗口加标题
（3）指定仿真时间
（4）设置示波器的显示刻度
解：（1）点击模型窗口的Format ，在下拉菜单中中有两个选项：Flip block 和Rotate block 。

Flip block 可使模块旋转180度（快捷键Ctrl+I ），Rotate block 可以模块顺时针旋转90度（快捷键Ctrl+R ）。

（2）在模型窗口的上方双击鼠标左键，会出现的一个文本输入框，编辑输入标题。

（3）单击Simulation ，在下拉菜单中选择configuration parameters ，打开仿真环境参数对话框，在Solver 选项的start time 设置仿真起止时间。

（4）双击示波器Scope 调出显示屏幕，选择Axes properties 选项，改变Y-min 和Y-max 可改变显示刻度。

6.3用Simulink 建立如下控制系统的仿真模型，并对系统进行阶跃响应仿真模型：
22()48
G s s s =++ 解：
6.4在Simulink 环境下，设计一个PID 控制器，实现下面被控制对象的控制，并观察选择不同的PID 参数时对控制效果的影响：
220()212G s s s =++ 其中，系统输入信号分别选择阶跃信号和正弦信号。

解：当输入信号为阶跃信号时：
当输入信号为正弦信号时：。

第六章 习题及参考答案一、习题1、已知一个由下列差分方程表示的系统，x(n)、y(n)分别表示该系统的输入、输出信号：)1(21)()2(61)1(65)(-+=-+--n x n x n y n y n y （1）画出该系统的直接型结构； （2）画出该系统的级联型结构； （3）画出该系统的并联型结构。

2、已知某系统的系统函数为：)6.09.01)(5.01()9.21)(1()(211211------++-+-+=z z z z z z z H 请画出该系统的级联型结构。

3、已知FIR 滤波器的单位脉冲响应为)(8.0)(5n R n h n =， （1）求该滤波器的系统函数； （2）画出该滤波器的直接型结构。

4、已知滤波器的系统函数为：3213218.09.09.018.04.16.01)(-------+-+--=zz z z z z z H 请画出该滤波器的直接型结构。

5、已知滤波器的系统函数为：)8.027.11)(5.01()44.11)(1(3)(211211------+--+--=z z z z z z z H 请画出该滤波器的级联型结构和并联型结构。

6、已知某因果系统的信号流图如下图所示：x(n)y(n)-25-3求该系统的系统函数和单位脉冲响应。

7、已知某系统的信号流图如下图所示：x(n)y(n)求该系统的系统函数和极点。

8、已知IIR 滤波器的系统函数为：4.035.04.046.16.14)(2323++++--=z z z z z z z H （1）画出级联型网络结构，要求利用MATLAB 分解H(z); （2）用MATLAB 验证所求的级联型结构是否正确。

9、已知IIR 滤波器的系统函数为：3213214.035.04.016.141.158.12.5)(-------++-++=zz z z z z z H （1）画出该系统的并联型网络结构，要求用MATLAB 分解； （2）用MATLAB 验证（1）中所求的并联型结构是否正确。

1．、2． 用冲激响应不变法将以下 )(s H a变换为变换为 )(z H ，抽样周期为T 。

为任意正整数 ,)()( )2()()()1(022n s s As H b a s as s H na a -=+++=分析：①冲激响应不变法满足)()()(nT h t h n h a nTt a===，T 为抽样间隔。

这种变换法必须)(s H a 先用部分分式展开。

②第（②第（22）小题要复习拉普拉斯变换公式1!][+=n n S n t L ,na n t s a S S As H t u n t Ae t h )()()()!1()(010-=⇔-=-，可求出可求出 )()()(kT Th t Th k h a kTt a===，|又 dz z dX z k kx )()(-⇔，则可递推求解。

解: (1)22111()()2a s a H s s a b s a jb s a jb ⎡⎤+==+⎢⎥+++++-⎣⎦[])( 21)()()(t u e e t h tjb a t jb a a --+-+= 由冲激响应不变法可得：由冲激响应不变法可得：由冲激响应不变法可得：[]()()()() ()2a jb nT a jb nT a T h n Th nT e e u n -+--==+ 11011() () 211naT jbT aT jbT n T H z h n z e e z e e z ∞------=⎡⎤==+⎢⎥--⎣⎦∑ 2211cos 21cos 1 ------+--⋅=z e bT z e bT z e T aT aT aT(2) 先引用拉氏变换的结论[]1!+=n n s n t L 可得：可得： n a s s As H )()(0-=))()!1()(10t u n t Ae t h n t s a -=-则 )()!1()()()(10k u n kT Ae T Tk Th k h n kT s a -⋅==- dzz dX zk kx az k u a ZZk )()( , 11)( 1-−→←-−→←-且按)11()()!1( )()!1( )()(111111000--∞=---∞=----=-==∑∑ze dz d z n AT e z k n T TA z k h z H T s n n kkT s n n k k可得⎪⎪⎩⎪⎪⎨⎧=-=-=•••---,3,2)1(1,1)(111000n z e z e AT n z e ATz H n T s T S n T s ，可以递推求得：2. 已知模拟二阶巴特沃思低通滤波器的归一化系统函数为：2'4142136.111)(s s s H a ++=而而3dB 截止频率为50Hz 的模拟滤波器，需将归一化的)('s H a 中的s 变量用502⨯πs 来代替424'108696044.928830.444108696044.9)100()(⨯++⨯==s s s H s H a a π：设系统抽样频率为设系统抽样频率为Hz f s 500=，要求从这一低通模拟滤波器设计一个低通数字滤波器，采用阶跃响应不变法。

第六章 习题解答(部分)［1］解：对采样数字系统，数字频率ω与模拟角频率Ω之间满足线性关系T Ω=ω。

因此，当时，ms T 01.0=TT cc 8πω==Ω，Hz T f c c 6251612==Ω=π 当s T µ5=时， TT c c 8πω==Ω，Hz T f c c 125001612==Ω=π［2］解：的极点为：，)(s H a jb a s +−=1jb a s −−=1将部分分式展开： )(s H a )(21)(21)(jb a s j jb a s js H a +−−−+−−−=所以有1)(1)(121121)(−+−−−−−−+−=z e j z e j z H T jb a T jb a通分并化简整理得：TT T e z bT e z bTe z z H ααα2211cos 21sin )(−−−−−−+−=［3］解：归一化原型低通滤波器与带通滤波器之间的频率变换关系为：B⋅ΩΩ−Ω=Ω22s rad p p /1002210×=ΩΩ=Ωπ，s rad B /2002×=π，dB p 2=δs rad s /80021×=Ωπ，s rad s /124022×=Ωπ，dB s 15=δ因此，归一化原型低通滤波器的通带频率p Ω取1，通带处最小衰减为2dB 。

同理可得归一化原型低通滤波器的阻带频率分别为：9375.31221=ΩΩ−Ω=ΩΩ=Ωs Bs ， 1597.62222=ΩΩ−Ω=ΩΩ=Ωs Bs因此，归一化原型低通滤波器的阻带频率9375.3),min(21=ΩΩ=Ωs s s ，这是因为取较小的频率值，则较大的频率处一定满足衰减要求，阻带处最大衰减为15dB 。

利用巴特沃斯低通滤波器设计归一化原型低通滤波器)(s H 利用归一化原型低通滤波器的指标，得巴特沃斯低通滤波器阶数N444.19372.31lg 2110110lg 5.12.0=⎟⎠⎞⎜⎝⎛⎟⎟⎠⎞⎜⎜⎝⎛−−≥N 取，查表的归一化巴特沃斯原型低通滤波器的系统函数 2=N 14142.11)(2++=s s s H LP由归一化原型低通滤波器变换到实际模拟带通滤波器22202220222)(4142.1)()()(202B s sB s s B s s H s H Bs s s LP BP +Ω++Ω+==⋅Ω+= ［4］解：（1）用冲激响应不变法① 确定数字滤波器指标rad p 3/πω=，dB p 3=δ rad s 5/4πω=，dB s 15=δ② 将数字滤波器指标转换为相应的模拟滤波器指标。