第5章 MATLAB在信号处理中的应用习题
MATLAB语言基础与应用(第二版)第5章 习题答案
第5章习题与答案5.1用矩阵三角分解方法解方程组123123123214453186920x x x x x x x x x +-=⎧⎪-+=⎨⎪+-=⎩ 解答：>>A=[2 1 -1;4 -1 3;6 9 -1] A =2 1 -1 4 -13 6 9 -1 >>b=[14 18 20]; b =14 18 20 >> [L, U, P]=lu(A) L =1.0000 0 0 0.6667 1.0000 0 0.3333 0.2857 1.0000 U =6.0000 9.0000 -1.0000 0 -7.0000 3.6667 0 0 -1.7143 P =0 0 1 0 1 0 1 0 0 >> y=backsub(L,P*b’) y =20.0000 4.6667 6.0000 >> x=backsub(U,y) x =6.5000 -2.5000 -3.5000 5.2 Cholesky 分解方法解方程组123121332352233127x x x x x x x ++=⎧⎪+=⎨⎪+=⎩ 解答：>> A=[3 2 3;2 2 0;3 0 12] A =3 2 32 2 03 0 12>> b=[5;3;7]b =537>> L=chol(A)L =1.7321 1.1547 1.73210 0.8165 -2.44950 0 1.7321>> y=backsub(L,b)y =-11.6871 15.7986 4.0415>> x=backsub(L',y)x =-6.7475 28.8917 49.93995.3解答:观察数据点图形>> x=0:0.5:2.5x =0 0.5000 1.0000 1.5000 2.0000 2.5000 >> y=[2.0 1.1 0.9 0.6 0.4 0.3]y =2.0000 1.1000 0.9000 0.6000 0.4000 0.3000 >> plot(x,y)图5.1 离散点分布示意图从图5.1观察数据点分布，用二次曲线拟合。
数字信号处理,第5章课后习题答案
第五章习题与上机题5.1 已知序列12()(),0 1 , ()()()nx n a u n a x n u n u n N =<<=--，分别求它们的自相关函数，并证明二者都是偶对称的实序列。
解：111()()()()()nn mx n n r m x n x n m a u n au n m ∞∞-=-∞=-∞=-=-∑∑当0m ≥时，122()1mmnx n ma r m aaa∞-===-∑ 当0m <时，122()1m mnx n a r m aaa -∞-===-∑ 所以，12()1mx ar m a =-2 ()()()()N x n u n u n N R n =--=22210121()()()()()1,0 =1,00, =()(1)x NN n n N mn N n m N r m x n x n m Rn R n m N m N m N m m Nm N m R m N ∞∞=-∞=-∞--=-=-=-=-⎧=--<<⎪⎪⎪⎪=-≤<⎨⎪⎪⎪⎪⎩-+-∑∑∑∑其他从1()x r m 和2()x r m 的表达式可以看出二者都是偶对称的实序列。
5.2 设()e()nTx n u n -=，T 为采样间隔。
求()x n 的自相关函数()x r m 。
解：解：()()()()e()e ()nTn m T x n n r m x n x n m u n u n m ∞∞---=-∞=-∞=-=-∑∑用5.1题计算1()x r m 的相同方法可得2e()1e m Tx Tr m --=-5.3 已知12()sin(2)sin(2)s s x n A f nT B f nT ππ=+，其中12,,,A B f f 均为常数。
求()x n 的自相关函数()x r m 。
解：解：()x n 可表为)()()(n v n u n x +=的形式，其中)2sin()(11s nT f A n u π=，=)(n v 22sin(2)s A f nT π，)(),(n v n u 的周期分别为 s T f N 111=，sT f N 221=，()x n 的周期N 则是21,N N 的最小公倍数。
MATLAB信号处理例题
◆例1设方波的数学模型为：]5sin 513sin 31[sin 4)(000 t t t E t f T ，基频： T20用MATLAB 软件完成该方波的合成设计◆ MATLAB 源程序t=-10:0.1:10; %设定一个数组有201个点，方波周期为20e=5;w=pi/10; %设定方波幅值为5，w 代表w0m=-5*sign(t); %给定幅值为5的方波函数y1=(-4*e/pi)*sin(w*t); %计算1次谐波y3=(-4*e/pi)*(sin(w*t)+sin(3*w*t)/3); %计算3次谐波y5=(-4*e/pi)*(sin(w*t)+sin(3*w*t)/3+sin(5*w*t)/5); %计算5次谐波plot(t,y1,'y');hold; grid; %用黄色点线画出1次谐波及网格线,并在同一张图上画其余曲线plot(t,y3,'g'); %用绿色点线画出3次谐波plot(t,y5,'b'); %用蓝色点线画出5次谐波plot(t,m,'-k'); %用黑色实线画方波title('方波合成');xlabel('t');ylabel('f(t)'); %为图形加上标题n=50; %合成任意次方波，n 决定方波的合成次数，在此给定50yn=0; %设置初始值for i=1:nyn=yn+(-4*e/pi)*(1/(2*i-1))*sin((2*i-1)*w*t);end; %计算n 次谐波合成plot(t,yn,'r') %用红色实线画出n 次谐波合成◆ 从图中我们可以看到Gibbs 现象。
在函数的间断点附近，增加傅里叶级数的展开次数，虽然可以使其间断点附近的微小振动的周期变小，但振幅却不能变小。
此现象在控制系统表现为：当求控制系统对阶跃函数的响应时，超调量总是存在的。
Matlab 在信号处理中的应用
Matlab 在信号处理中的应用1.产生锯齿波或三角波信号函数sawtooth格式一：x=sawtooth（t）功能：产生周期为2pi，振幅从-1到1的锯齿波。
在2pi的整数倍处值为-1到1，这一段波形斜率为1/pi。
格式二：sawtooth（t，width）功能：产生三角波，width在0到1 之间。
例：产生周期为0.02的三角波。
Fs=10000；t=0:1/Fs：1X1=sawtooth（2*pi*50*t，0）；X2=sawtooth（2*pi*50*t，1）；Subplot（2,1,1），plot（t，x1），axis（[0,0.2，-1,1]）；Subplot（2,1,2），plot（t，x2），axis（[0,0.2，-1,1]）；2.产生方波信号函数square格式：x=square（t）功能：产生周期为2pi，振幅从-1到1的方波。
格式：x=square（t，duty）功能：产生指定周期的方波，duty为正半周期的比例。
例：产生周期为0.02的方波。
Fs=10000；t=0:1/Fs：1X1=square（2*pi*50*t，20）；X2=square（2*pi*50*t，80）；Subplot（2,1,1），plot（t，x1），axis（[0,0.2，-1.5,1.5]）；Subplot（2,1,2），plot（t，x2），axis（[0,0.2，-1.5,1.5]）；3.产生sinc函数波形函数sinc格式：y=sinc（x）功能：sinc（x）用于计算sinc函数，即：1 t=0Sinc（t）= sin（pi*t）/（pi*t）Sinc函数十分重要，它的傅里叶变换正好是振幅为1的矩阵脉冲例：产生sinc函数波形。
X=linspace（-4,4）Y=sinc（x）Plot（x，y）4.产生非周期方波信号函数rectpuls格式：y=rectpuls（t）功能：产生非周期方波信号，方波的宽度为时间轴的一半格式：y=rectpuls（t，w）功能：产生指定宽度为w的非周期方波。
Matlab在信号处理中的应用
MATAB在信号处理中的应用绪言：MATLAB是由美国mathworks公司发布的主要面对科学计算、可视化以及交互式程序设计的高科技计算环境。
它将数值分析、矩阵计算、科学数据可视化以及非线性动态系统的建模和仿真等诸多强大功能集成在一个易于使用的视窗环境中，为科学研究、工程设计以及必须进行有效数值计算的众多科学领域提供了一种全面的解决方案，并在很大程度上摆脱了传统非交互式程序设计语言（如C、Fortran）的编辑模式，代表了当今国际科学计算软件的先进水平。
MATLAB和Mathematica、Maple、MathCAD并称为四大数学软件。
它在数学类科技应用软件中在数值计算方面首屈一指。
MATLAB可以进行矩阵运算、绘制函数和数据、实现算法、创建用户界面、连接其他编程语言的程序等，主要应用于工程计算、控制设计、信号处理与通讯、图像处理、信号检测、金融建模设计与分析等领域。
MATLAB的基本数据单位是矩阵，它的指令表达式与数学、工程中常用的形式十分相似，故用MATLAB来解算问题要比用C，FORTRAN等语言完成相同的事情简捷得多，并且MATLAB也吸收了像Maple等软件的优点，使MATLAB成为一个强大的数学软件。
在新的版本中也加入了对C，FORTRAN，C++，JAVA的支持。
可以直接调用,用户也可以将自己编写的实用程序导入到MATLAB函数库中方便自己以后调用，此外许多的MATLAB爱好者都编写了一些经典的程序，用户可以直接进行下载就可以用。
MATLAB是一个包含大量计算算法的集合。
其拥有600多个工程中要用到的数学运算函数，可以方便的实现用户所需的各种计算功能。
函数中所使用的算法都是科研和工程计算中的最新研究成果，而前经过了各种优化和容错处理。
在通常情况下，可以用它来代替底层编程语言，如C和C++ 。
在计算要求相同的情况下，使用MATLAB的编程工作量会大大减少。
MATLAB的这些函数集包括从最简单最基本的函数到诸如矩阵，特征向量、快速傅立叶变换的复杂函数。
数字信号处理课后习题Matlab作业
数字信号处理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程序，计算下列序列的循环卷积。
实验五 在数字信号处理中Matlab的应用(答案)
1、程序：>> p=8;q=2; t=0:16; x=exp(-(t-p).^2./q); subplot(231) stem(t,x); xlabel('n');ylabel('Xa');title('p=8,q=2'); hold on N=16;fs=1; n=0:N-1;t=n/fs; x=exp(-(t-p).^2./q); y=fft(x,N); mag=abs(y); hou=mag.*2./N; f=n*fs/N; subplot(234) plot(f,hou); xlabel('频率/Hz');ylabel('幅值');title('p=8,q=2'); >> p=8;q=4; t=0:16; x=exp(-(t-p).^2./q); subplot(232) stem(t,x); xlabel('n');ylabel('Xa');title('p=8,q=4'); hold on N=16;fs=1; n=0:N-1;t=n/fs; x=exp(-(t-p).^2./q); y=fft(x,N); mag=abs(y); hou=mag.*2./N; f=n*fs/N; subplot(235) plot(f,hou); xlabel('频率/Hz');ylabel('幅值');title('p=8,q=4'); >> p=8;q=8; t=0:16; x=exp(-(t-p).^2./q); subplot(233) stem(t,x); xlabel('n');ylabel('Xa');title('p=8,q=8'); hold on N=16;fs=1; n=0:N-1;t=n/fs;
Matlab在信号处理中的应用63858-
· function [返回变量列表] 函数名(输入变量列表) 注释说明语句段, 由% 引导 输入、返回变量格式的检测 函数体语句
2 MATLAB在采样与波形发生的应用
· 数字信号处理的对象，是在采样时钟的控制之下，通过A/D 转换器一定的采样频率对模拟信号进行采样得到的。根据采 样定理，采样率必须大于模拟信号的最高采样率（Nyquist ）的2倍。但是在许多情况下，要求信号以不同的频率采样 ，改变采样时钟虽然可行，但是并不可取。这时需要对采样 数据进行处理：或者用抽取的方法降低其采样率（下采样） ；或者用内插的方法提高其采样率（上采样）；或者两者兼 有之（重采样）。
· for 循环变量=s1:s3:s2 循环体语句组
end · 注意，这里的循环语句是以end 结尾的，这和C 语言的结构不完全一致。 · 试探结构： · MATLAB 从 5.2 版本开始提供了一种新的试探式语句结构，其一般的形式为： · try
语句段1 catch
语句段2 end · 本语句结构首先试探性地执行语句段 1,如果在此段语句执行过程中出现错误，则 将错误信息赋给保留的lasterr 变量，并放弃这段语句，转而执行语句段 2 中的 语句。这种新的语句结构是C 等语言中所没有的。
1.5 MATLAB 的语句流程与控制
· 开关结构： · MATLAB 从 5.0 版开始提供了开关语句结构，其基本语句结构为：
· switch 开关表达式
case 表达式1 语句段1
case {表达式2,表达式3,..., 表达式m} 语句段2
... otherwise
语句段n end · MATLAB 开关语句与C 有区别： · 当开关表达式的值等于表达式 1 时，将执行语句段 1,执行完语句段 1 后将转出开关体 ，无需像C 语言那样在下一个case 语句前加break 语句，所以本结构在这点上和C 语 言是不同的。
数字信号处理第5章答案
第５章 无限脉冲响应(IIR)数字滤波器的设
计
特别是以理想滤波器特性作为Hd(ejω)时， 为了使ε2最小，
优化过程尽可能逼近Hd(ejω)的间断特性（即使过渡带最窄）， 而使通带出现较大过冲、 阻带最小衰减过小， 不能满足工
H(ejω)=|H(ejω)|ejθ(ω) 其中， |H(ejω)|称为幅频特性函数， θ(ω)称为相频特性函数。
常用的典型滤波器|H(ejω)|是归一化的， 即|H(ejω)|max=1， 下 的讨论一般就是针对归一化情况的。 对IIR数字滤波器， 通
常用幅频响应函数|H(ejω)|来描述设计指标， 而对线性相位特 性的滤波器， 一般用FIR数字滤波器设计实现。
计
图5.1.6
第５章 无限脉冲响应(IIR)数字滤波器的设
计
5.1.4 IIR-DF的直接设计法
所谓直接设计法， 就是直接在数字域设计IIR[CD*2]DF 的方法。 相对而言， 因为从AF入手设计DF是先设计相应的 AF， 然后再通过s-z平面映射， 将Ha(s)转换成H(z)， 所以 这属于间接设计法。 该设计法只能设计与几种典型AF相对 应的幅频特性的DF。 而需要设计任意形状幅频特性的DF时， 只能用直接设计法。 直接设计法一般都要借助于计算机进行 设计， 即计算机辅助设计（CAD）。 现在已有多种DF优化 设计程序。 优化准则不同， 所设计的滤波器特点亦不同。所 以最主要的是建立优化设计的概念， 了解各种优化准则的 特点， 并根据设计要求， 选择合适的优化程序设计DF。
≤≤
(5.1.1)
≤
(5.1.2)
MATLAB信号处理例题
◆例1设方波的数学模型为：]5sin 513sin 31[sin 4)(000 +++=t t t E t f T ωωωπ ，基频： Tπω20=用MATLAB 软件完成该方波的合成设计◆ MATLAB 源程序t=-10:0.1:10; %设定一个数组有201个点，方波周期为20e=5;w=pi/10; %设定方波幅值为5，w 代表w0m=-5*sign(t); %给定幅值为5的方波函数y1=(-4*e/pi)*sin(w*t); %计算1次谐波y3=(-4*e/pi)*(sin(w*t)+sin(3*w*t)/3); %计算3次谐波y5=(-4*e/pi)*(sin(w*t)+sin(3*w*t)/3+sin(5*w*t)/5); %计算5次谐波 plot(t,y1,'y');hold; grid; %用黄色点线画出1次谐波及网格线,并在同一张图上画其余曲线plot(t,y3,'g'); %用绿色点线画出3次谐波plot(t,y5,'b'); %用蓝色点线画出5次谐波plot(t,m,'-k'); %用黑色实线画方波title('方波合成');xlabel('t');ylabel('f(t)'); %为图形加上标题 n=50; %合成任意次方波，n 决定方波的合成次数，在此给定50yn=0; %设置初始值for i=1:nyn=yn+(-4*e/pi)*(1/(2*i-1))*sin((2*i-1)*w*t);end; %计算n 次谐波合成plot(t,yn,'r') %用红色实线画出n 次谐波合成◆ 从图中我们可以看到Gibbs 现象。
在函数的间断点附近，增加傅里叶级数的展开次数，虽然可以使其间断点附近的微小振动的周期变小，但振幅却不能变小。
此现象在控制系统表现为：当求控制系统对阶跃函数的响应时，超调量总是存在的。
MATLAB在信号处理中的应用
基本信号的MATLAB实现：1.δ函数：x=[zeros(1,10),1,zeros(1,10)];stem(x);%注：若是产生序列可以用stem(x)代替stem(n,x);当然通过改变为值为1时的序列下标，实现δ函数的时移；2.ε函数：n=-20:20;t=(n>=0);stem(t);%通过改变(n>=i)来实现ε函数的时移；3.单位斜坡函数：n=-20:20;t=n.*(n>0);stem(t);%通过改变n.*(n>i)或(n-i).*(n>i)来实现不同的单位斜坡函数时移；4.复指数序列：>> n=-10:10;t=0.1+j*0.3;x=exp(t*n);>> subplot(221);stem(real(x));>> subplot(222);stem(imag(x));>> subplot(223);stem(abs(x));>> subplot(224);stem((180/pi*angle(x)));5.随机序列：使用rand(1,n)和randn(1,n)产生随机序列；6.MATLAB信号工具箱还提供了一些其他的常用信号，如：SQUARE,SAWTOOTH,SINC,DIRIC,DIRICHLET,RECTPULS和PULSTRAN，具体用法参考help文件。
MATLAB常用函数：real(x):返回复数的实部；imag(x):返回复数的虚部；abs(x):返回复数的模；angle(x):返回复数的相角；rand(1,n):返回长度为n的[0,1]上均匀分布的随机序列；randn(1,n):返回长度为n的均值为0，方差为1的高斯随机序列，即白噪声序列；信号的基本运算：1.信号相加：function [y,n]=sig_add(x1,n1,x2,n2)% Implements y(n)=x1(n)+x2(n);%x1,x2:序列%n1,n2:序列的起始/终止下标n=min(min(n1),min(n2)):max(max(n1),max(n2));y1=zeros(1,length(n));y2=y1;y1(find((n>=min(n1))&(n<=max(n1))==1))=x1;y2(find((n>=min(n2))&(n<=max(n2))==1))=x2;y=y1+y2;2.信号相乘：function [y,n]=sig_mult(x1,n1,x2,n2)% Implements y(n)=x1(n)*x2(n);%x1,x2:序列%n1,n2:序列的起始/终止下标n=min(min(n1),min(n2)):max(max(n1),max(n2));y1=zeros(1,length(n));y2=y1;y1(find((n>=min(n1))&(n<=max(n1))==1))=x1;y2(find((n>=min(n2))&(n<=max(n2))==1))=x2;y=y1.*y2;3.信号移位：function [y,n]=sig_shift(x,m,n0)% Implements y(m+k)=x(n)n=m+n0;y=x;4.序列折叠：function [y,n]=sig_fold(x,m)% Implements y(n)=x(-n)y=fliplr(x);n=-fliplr(n)5.序列奇偶性：function[xeven,xodd,m]=sig_even_odd(x,n)%Real signal decomposition into even and odd parts if(imag(x)~=0)error('x is not a real sequence');endm=-fliplr(n);m1=min([m,n]);m2=max([m,n]);m=m1:m2;nm=n(1)-m(1);n1=1:length(n);x1=zeros(1,length(m));x1(n1+nm)=x;x=x1;xeven=0.5*(x+fliplr(x));xodd=0.5*(x-fliplr(x));下面的用于求奇偶分量的函数的特点在于：它根据输入的序列的中点为对称中心而求奇偶分量：function [xeven,xodd]=circevod(x)if any(imag(x)~=0)error('The sequence is not real.')endN=length(x);n=0:(N-1);xeven=0.5*(x+x(sigmod(-n,N)+1));xodd=0.5*(x-x(sigmod(-n,N)+1));例：n=-5:17;>> x=10*0.8.^n;;>> [xeven,xodd]=circevod(x);>> subplot(311),stem(x); >> subplot(312),stem(xeven); >> subplot(313),stem(xodd);020400510152025020400510152025-200206. 序列圆周移位： function y=cirshift(x,m,N) %x:input sequence %m:shift number %N:show length if length(x)>Nerror('N must be greater then length(x)'); endx=[x zeros(1,N-length(x))];n=[0:N-1];n=sigmod(n-m,N);y=x(n+1);sigmod函数用来找出周期序列任意位置n所对应的主值有限序列x(n)中的位置m(若周期序列由有限长序列x(n)产生，周期为N)function m=sigmod(n,N)m=rem(n,N);m=m+N;m=rem(m,N);例：由下例可以看出，循环移位不会改变频谱的幅值，只会改变其相位。
第5章MATLAB在信号处理中的应用
5.2.3 序列翻褶与序列累加运算
序列翻褶：y(n)=x(-n)。MATLAB可实现: y=fliplr(x) 序列累加的数学描述为：y ( n ) MATLAB实现：y=cumsum(x)
i ns
x (i )
n
第 5 章 MATLAB在信号处理中的应用
5.2.4 两序列的卷积运算
两序列卷积运算： y(n) x1 (n) x2 (n)
3．零－极点增益模型 连续系统： 离散系统：
H ( z) k
( z q1 )( z q2 ) ( z qM ) ( z p1 )( z p2 ) ( z pN )
第 5 章 MATLAB在信号处理中的应用
4．极点留数模型 连续系统： 离散系统：
H ( s) r r1 r 2 N s p1 s p2 s pN
第 5 章 MATLAB在信号处理中的应用
5.4 线性时不变系统
5.4.1 系统的描述
1．常系数线性微分/差分方程
y
(N)
(t ) ai y (t ) bi x(i ) (t )
(i ) i 0 i 0
N 1
M
2．系统传递函数 连续系统： 离散系统：
Y (s) bM s M bM 1s M 1 b1s b0 H ( s) N X ( s) s aN 1s N 1 a1s a0
第 5 章 MATLAB在信号处理中的应用
3．实指数序列
x(n) a n
n
aR
直接实现：n=[ns:nf]; x=a.^n; 函数实现：利用实指数序列 的生成函数rexpseq，即 function [x, n] = rexpseq(a, ns, nf) n = [ns:nf]; x = a.^n; 4．复指数序列
Matlab在信号处理中的应用范例
1.2 MATLAB特色举例
· >> [x,y] = meshgrid(-3:1/8:3); z = 3*(1-x).^2.*exp(-(x.^2) - (y+1).^2)- 10*(x/5 - x.^3 -
y.^5)... .*exp(-x.^2-y.^2)- 1/3*exp(-(x+1).^2 - y.^2); surf(x,y,z), shading interp; colorbar
X1
1.6 MATLAB 函数编写方法与应用
· 所 谓 MATLAB 程序，大致分为两类：M 脚本文件(M-script) 和 M 函数(M-function), 它们均是普通的ASCII 码构成的文件。 M 脚本文件中包含一族由MATLAB 语言所支持的语句，它类似 于 DOS 下的批处理文件，它的执行方式很简单，用户只需在 MATLAB 的提示符>> 下键入该M 文件的文件名，这样 MATLAB 就会自动执行该M 文件中的各条语句，并将结果直接返回到 MATLAB 的工作空间。M 函数格式是MATLAB 程序设 计的主流，一般情况下， 不建议您使用M 脚本文件格式编程。 · MATLAB 的M 函数是由function 语句引导的，其基本格式 如下： · function [返回变量列表] = 函数名(输入变量列表) 注释说明语句段, 由% 引导 输入、返回变量格式的检测 函数体语句
MATLAB 还支持其他运算，如取整、求余数等。可以使用rond（）, fix(), rem() 等来实现。
X1
1.5 MATLAB 的语句流程与控制
· 作为一种常用的编程语言，MATLAB 支持各种流程控制结构，如循环结构、条件 转移结构、客观结构等另外MATLAB 还支持一种新的结构--- 试探结构。 · 循环语句有两种结构： for ... end 结构和while ... end 结构。 这两种语句结构不 完全相同，各有各的特色。for ... end 语句通常的调用格式为： ·for 循环变量=s1:s3:s2循
MATLAB在信号与系统中的应用(习题)
实验四MATLAB在信号与系统中的应用第一部分连续时间信号与系统（习题）2、编写程序，完成下列连续信号波形1）F(t)=u(t-2)+u(t-4) (-2≤t≤6)2）F(t)=e-t cos(4πt) (0≤t≤3)3）F(t)=3e(0.2+j0.5π)t (0≤t≤4π)1）方法一>> cleart0=-2;tf=6;dt=0.05;t1=2;t2=4;t=[t0:dt:tf];st=length(t);n1=floor((t1-t0)/dt);x1=[zeros(1,n1),ones(1,st-n1)];n2=floor((t2-t0)/dt);x2=[zeros(1,n2),ones(1,st-n2)];x=x1+x2;stairs(t,x),grid onaxis([0,5,0,3])方法二>> cleart0=-2;tf=6;dt=0.05;t1=2;t2=4;t=[t0:dt:tf];st=length(t);x1=(t>=t1);x2=(t>=t2);x=x1+x2;stairs(t,x),grid on2）cleart0=0;tf=3;dt=0.05;t1=4*pi;t=[t0:dt:tf];st=length(t);alpha=-1;w=0;x1=exp((alpha+j*w).*t);x2=cos(t1.*t);x=x1.*x2;subplot(1,1,1),plot(t,x),grid on3）cleart0=0;tf=4*pi;dt=0.05;t=[t0:dt:tf];st=length(t);alpha=0.2;w=0.5*pi;x=3*exp((alpha+j*w).*t);subplot(2,1,1),plot(t,real(x)),gridon, title ' real x'subplot(2,1,2),plot(t,imag(x)),gridon, title ' imag x'3、已知某连续时间系统的微分方程为：()2()()()2()r t r t r t e t e t ''''++=+，求当输入信号为2()()t e t e t ε-=，该系统的零状态响应r(t)。
第5章MATLAB在信号处理中的应用
信号的MATLAB表示及运算
本节首先讨论用MATLAB表示和分析连 续信号的问题 严格说来，只有用符号推理的方法才能分 析连续系统。用数值方法是不能表示连续 信号的，因为它给出的是各个样本点的数 据。只有当样本点取得很密时才可看成连 续信号。所谓密，是相对于信号变化的快 慢而言。以下均假定相对于采样点密度而 言，信号变化足够慢
plot(t,ft)
axis([0,4,-0.5,1.5])
(2)方法二： 自编产生以指定位置为中心,指定宽度的矩 形脉冲
clc,clear t=-5:0.001:5; delta=[t>=-1&t<=1]; %以0为中心，门宽为2的矩形脉冲 plot(t,delta) axis([-5,5,-0.1,1.1]) grid
在MATLAB中能表示为 k=[-3,-2,-1,0,1,2,3,4]; f=[2,1,-1,0,1,4,3,7];
stem(k,f,’filled’)
说明：位置向量k是顺序增加的整数，所以在 MATLAB中可以简化表示为k=-3:4；有些序列，因为 采样位置从k=0开始，就省写了k向量，只用x向量来表 示序列x[k]
例如，产生一个幅度为1，基频为4Hz，占空比 为20%的周期方波
A=1;f=4; w=2*pi*f;T=1/f t=0:0.001:10*T; y=A*square(w*t,50);
plot(t,y);
xlabel('Time(s)'); title('square wave') axis([0,2.6,-1.5,1.5])
二.常用的典型序列
单位脉冲序列 单位阶跃序列 复指数序列
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
1、根据下图实际电路知识，可求解其方程为：y ’’(t)+3y ’(t)+6y(t)=6x(t)，要求用impulse 函数求解单位冲激响应、用step 函数求单位阶跃响应。
仿真程序：ts = 0;te = 5;dt = 0.01; sys = tf([6],[1,3,6]); t = ts:dt:te;y = impulse(sys,t); plot(t,y);xlabel('Time(sec)') ylabel('h(t)') untiled21.mts = 0;te = 4;dt = 0.01; sys = tf([6],[1,3,6]); t = ts:dt:te; y = step(sys,t); plot(t,y);xlabel('Time(sec)') ylabel('g(t)')2、 已知一个因果系统的差分方程为()()()()()()11111212336226y n y n x n x n x n x n +-=+-+-+-满足初始条件()10,x (1)0,y -=-=求系统的单位样值响应和阶跃响应。
编写MATLAB 程序如下：a=[1,0,1/3,0]; %x(n)系数b=[1/6,1/2,1/2,1/6]; %y(n)系数 N=32; %取32点作图 n=0:N-1;hn=impz(b,a,n); %求单位样值响应 gn=dstep(b,a,n); %求阶跃响应 subplot(1,2,1),stem(n,hn,'k'); title(‘系统的单位样值响应');ylabel('h(n)');xlabel('n');axis([0,N,-1.1*min(hn),1.1*max(hn)]); subplot(1,2,2),stem(n,gn,'k'); title('系统的阶跃响应'); ylabel('g(n)');xlabel('n');axis([0,N,-1.1*min(gn),1.1*max(gn)]);3、 已知一个因果系统的差分方程为6()2(4)()3(2)3(4)(6)y n y n x n x n x n x n --=--+---满足初始条件(1)0,(1)0,y x -=-=求系统单位冲激响应和单位阶跃响应。
时间轴上N 取32点作图。
编写MATLAB 程序如下：x01=0;y01=0;N=32; %赋初始条件和采样点数 a=[1,0,0,0,-1/3,0,0]; %输入差分方程系数 b=[1/6,0,-1/2,0,1/2,0,-1/6];xi=filtic(b,a,0); %求等效初始条件的输入序列 n=0:N-1; %建立N 点的时间序列x1=[n==0]; %建立输入单位冲激信号进行滤波，求冲激响应 hn=filter(b,a,x1,xi); %对输入单位冲激信号就 x2=[n>=0];gn=filter(b,a,x2,xi);subplot(1,2,1),stem(n,hn); title('系统单位冲激响应'); subplot(1,2,2),stem(n,gn); title('系统单位阶跃响应');4、已知一个连续时间信号, 001()sin (2)sin (6)3f t f t f t ππ=+,01f H z =，取最高有限带宽频率05m f f =,分别显示原连续时间信号波形和2s m F f >、2s m F f =、2s m F f <三种情况下抽样信号的波形。
MATLAB 程序如下：dt=0.1;f0=1;T0=1/f0; fm=5*f0;Tm=1/fm; t=-2:dt:2;f=sin(2*pi*f0*t)+1/3*sin(6*pi*f0*t); subplot(4,1,1),plot(t,f);axis([min(t) max(t) 1.1*min(f) 1.1*max(f)]); title('原连续信号和抽样信号'); for i=1:3;fs=i*fm;Ts=1/fs; n=-2:Ts:2;f=sin(2*pi*f0*n)+1/3*sin(6*pi*f0*n); subplot(4,1,i+1),stem(n,f,'filled');axis([min(n) max(n) 1.1*min(f) 1.1*max(f)]); end5、求解题5中原连续信号波形和2s m F f >、2s m F f =、2s m F f <三种情况下的抽样信号波形所对应的幅度谱。
MATLAB 程序如下： dt=0.1;f0=1;T0=1/f0; t=-2:dt:2; N=length(t);f=sin(2*pi*f0*t)+1/3*sin(6*pi*f0*t); fm=5*f0;Tm=1/fm; wm=2*pi*fm; k=0:N-1; w1=k*wm/N;F1=f*exp(-j*t'*w1)*dt;subplot(4,1,1),plot(w1/(2*pi),abs(F1));axis([0 max(4*fm) 1.1*min(abs(F1)) 1.1*max(abs(F1))]); for i=1:3;if i<=2 c=0,else c=1,end fs=(i+c)*fm;Ts=1/fs; n=-2:Ts:2;f=sin(2*pi*f0*n)+1/3*sin(6*pi*f0*n); N=length(n); wm=2*pi*fs; k=0:N-1; w=k*wm/N;F=f*exp(-j*n'*w)*Ts;subplot(4,1,i+1),plot(w/(2*pi),abs(F));axis([0 max(4*fm) 1.1*min(abs(F)) 1.1*max(abs(F))]); end6、用时域卷积推导出的内插公式重建题4给定的信号。
001()sin (2)sin (6)3f t f t f t ππ=+MATLAB 程序如下： f0=1;T0=1/f0;dt=0.01; fm=5*f0;Tm=1/fm; t=0:dt:3*T0;x=sin(2*pi*f0*t)+1/3*sin(6*pi*f0*t);subplot(4,1,1),plot(t,x);axis([min(t) max(t) 1.1*min(x) 1.1*max(x)]); title('用时域卷积重建抽样信号'); for i=1:3;fs=i*fm;Ts=1/fs; n=0:(3*T0)/Ts t1=0:Ts:3*T0;x1=sin(2*pi*n*f0/fs)+1/3*sin(6*pi*n*f0/fs);T_N=ones(length(n),1)*t1-n'*Ts*ones(1,length(t1)); xa=x1*sinc(fs*pi*T_N); subplot(4,1,i+1),plot(t1,xa);axis([min(t1) max(t1) 1.1*min(xa) 1.1*max(xa)]); end7、求[]()0,1,2,3,4,5,6,7x n =，07n ≤≤的DTFT ，将（ 2π-，2π）区间分成500份。
要求：（1）画出原信号。
（2）画出由离散时间傅里叶变换求得的幅度谱()j X eω和相位谱a rg [()]j X eω图形。
MATLAB 程序如下：xn=[0,1,2,3,4,5,6,7]; N=length(xn); n=0:N-1;w=linspace(-2*pi,2*pi,500); X=xn*exp(-j*n'*w);subplot(3,1,1),stem(n,xn,'k'); ylabel('x(n)');subplot(3,1,2),plot(w,abs(X),'k')axis([-2*pi,2*pi,1.1*min(abs(X)),1.1*max(abs(X))]); ylabel('幅度谱');subplot(3,1,3),plot(w,angle(X),'k')axis([-2*pi,2*pi,1.1*min(angle(X)),1.1*max(angle(X))]); ylabel('相位谱');8、已知[]7,6,5,4,3,2,1,0)(=n x ，求)(n x 的DFT 和IDFT 。
要求： (1) 画出序列傅里叶变换对应的)(k X 和[])(arg k X 图形。
(2) 画出原信号与傅里叶逆变换IDFT [])(k X 图形进行比较。
解：MATLAB 程序如下： xn=[0,1,2,3,4,5,6,7]; N=length(xn);n=0:N-1;k=0:N-1;Xk=xn*exp(-j*2*pi/N).^(n'*k); x=(Xk*exp(j*2*pi/N).^(n'*k))/N; subplot(2,2,1),stem(n,xn);title('x(n)');subplot(2,2,2),stem(n,abs(x)); title('IDFT|X(k)|');subplot(2,2,3),stem(k,abs(Xk)); title('|X(k)|');subplot(2,2,4),stem(k,angle(Xk)); title('arg|X(k)|');9、采用脉冲响应不变法设计一个巴特沃斯数字低通滤波器，要求：πω25.0=p ，dBRp1=；πω4.0=s ，dBAs15=，滤波器采样频率2000=sF Hz。
解编写程序如下：%脉冲响应不变法设计巴特沃斯低通滤波器%数字滤波器指标wp=0.25*pi; %滤波器的通带截止频率ws=0.4*pi; %滤波器的阻带截止频率Rp=1;As=15; %输入滤波器的通阻带衰减指标ripple=10^(-Rp/20); %计算通带衰减对应的幅度值Attn=10^(-As/20); %计算阻带衰减对应的幅度值%转换为模拟滤波器指标Fs=2000;T=1/Fs; Omgp=wp*Fs;Omgs=ws*Fs; %模拟原型滤波器计算[n,Omgc]=buttord(Omgp,Omgs,Rp,As,'s')%计算阶数n和截止频率[z0,p0,k0]=buttap(n); %设计归一化的模拟原型滤波器ba1=k0*real(poly(z0)); %求原型滤波器系数baa1=real(poly(p0)); %求原型滤波器系数a[ba,aa]=lp2lp(ba1,aa1,Omgc);%变换为模拟低通滤波器%用脉冲响应不变法计算数字滤波器系数[bd,ad]=impinvar(ba,aa,Fs)[C,B,A]=dir2par(bd,ad) %转换成并联型%求数字系统的频率特性[H,w]=freqz(bd,ad);dbH=20*log10(abs(H)/max(abs(H))); %化为分贝值subplot(2,2,1),plot(w/pi,abs(H));ylabel('幅度');xlabel('频率/\pi');axis([0,1,0,1.1]);set(gca,'XTickMode','manual','XTick',[0,0.25,0.4,1]);set(gca,'YTickMode','manual','YTick',[0,Attn,ripple,1]);grid subplot(2,2,2),plot(w/pi,angle(H)/pi);ylabel('相位');xlabel('频率/\pi');axis([0,1,-1,1]);set(gca,'XTickMode','manual','XTick',[0,0.25,0.4,1]); set(gca,'YTickMode','manual','YTick',[-1,0,1]);grid subplot(2,2,3),plot(w/pi,dbH);ylabel('幅度(dB)');xlabel('频率/\pi');axis([0,1,-40,5]);set(gca,'XTickMode','manual','XTick',[0,0.25,0.4,1]); set(gca,'YTickMode','manual','YTick',[-50,-15,-1,0]);grid subplot(2,2,4),zplane(bd,ad);%axis([-1.1,1.1,-1.1,1.1]);ylabel('零极图');10、采用双线性变换法设计一个巴特沃斯数字低通滤波器，要求：πω25.0=p ，dBRp1=；πω4.0=s ，dBAs15=，滤波器的采样频率100=sF Hz。