MATLAB在数字信号处理中的应用.

1第7章MATLAB在数字信号处理中的应用27.1时域离散信号的产生及时域处理时域离散信号用x(n)表示，时间变量n（表示采样位置）只能取整数。因此，x(n)是一个离散序列，以后简称序列。用一个向量x不足以表示序列值x(n)。必须再用另一个等长的定位时间变量n。x和n同时使用才能完整地表示一个序列，由于n序列是按整数递增的，可简单地用其初值ns决定，因为它的终值nf取决于ns和x的长度length(x)，故可写成：n=[ns:nf]或n=[ns:nslength(x)1]3单位脉冲序列δ(n-n0)的生成函数impseqfunction[x,n]=impseq(n0,ns,nf)n=[ns:nf];x=[(n-n0)==0];单位阶跃序列u(n-n0)的生成函数stepseqfunction[x,n]=stepseq(n0,ns,nf)n=[ns:nf];x=[(n-n0)=0];P173表7.14例7.1序列的相加和相乘给出两个序列x1(n)和x2(n)。x1=[0,1,2,3,4,3,2,1,0];n1=[-2:6];x2=[2,2,0,0,0,-2,-2];n2=[2:8];要求它们的和ya及乘积yp。解：编程的思路是把序列长度延拓到覆盖n1和n2的范围，这样才能把两序列的时间变量对应起来，然后进行对应元素的运算。5%第七章例7.1q701离散信号序列的相加和相乘x1=[0,1,2,3,4,3,2,1,0];ns1=-2;%给定x1及ns1x2=[2,2,0,0,0,-2,-2];ns2=2;%给定x2及ns2nf1=ns1+length(x1)-1;nf2=ns2+length(x2)-1;ny=min(ns1,ns2):max(nf1,nf2)%y(n)的时间变量xa1=zeros(1,length(ny));xa2=xa1;%延拓序列初始化xa1(find((ny=ns1)&(ny=nf1)==1))=x1%给xa1赋值x1xa2(find((ny=ns2)&(ny=nf2)==1))=x2%给xa2赋值x2ya=xa1+xa2%序列相加yp=xa1.*xa2%序列相乘subplot(4,1,1),stem(ny,xa1,'.')%绘图subplot(4,1,2),stem(ny,xa2,'.')line([ny(1),ny(end)],[0,0])%画x轴subplot(4,1,3),stem(ny,ya,'.')line([ny(1),ny(end)],[0,0])%画x轴subplot(4,1,4),stem(ny,yp,'.')line([ny(1),ny(end)],[0,0])%画x轴6ny=-2-1012345678xa1=01234321000xa2=000022000-2-2ya=012365210-2-2yp=000086000007例7.2序列的合成和截取用例6.13的结果编写产生矩形序列RN(n)的程序。序列起点为n0，矩形序列起点为n1，长度为N（n0，n1，N由键盘输入）。并用它截取一个复正弦序列exp(jπn/8)，最后画出波形。解：建模:矩形序列可看成两个阶跃序列之差。即用MATLAB逻辑关系运算产生矩形序列x2(n)。而用它截取任何序列相当于元素群相乘x2.*x，也称为加窗运算。序列的合成和截取就是相加和相乘。1111()()()()NxnRnnUnnUnnN8%第七章例7.2程序q702离散信号序列的相加和相乘clear;closealln0=input('输入序列起点:n0=');%书上取-6N=input('输入序列长度:N=');%书上取15n1=input('输入位移:n1=');%书上取3n=n0:n1+N+5;%生成自变量数组u=[(n-n1)=0];%产生单位阶跃序列（u(n-n1)）x1=[(n-n1)=0]-[(n-n1-N)=0];%用阶跃序列差产生矩形序列x2=[(n=n1)&(n(N+n1))];%用逻辑式产生矩形序列x3=exp(j*n*pi/8).*x2;%对复正弦序列加矩形窗（元素群乘）subplot(2,2,1);stem(n,x1,'.');xlabel('n');ylabel('x1(n)');%标注axis([n0,max(n),0,1]);%定坐标范围subplot(2,2,3);stem(n,x2,'.');xlabel('n');ylabel('x2(n)');%标注axis([n0,max(n),0,1]);%定坐标范围9subplot(2,2,2);stem(n,real(x3),'.');xlabel('n');ylabel('x3(n)的实部');%标注line([n0,max(n)],[0,0]);%画横轴axis([n0,max(n),-1,1]);%定坐标范围subplot(2,2,4);stem(n,imag(x3),'.');xlabel('n');ylabel('x3(n)的虚部');%标注line([n0,max(n)],[0,0]);%画横轴axis([n0,max(n),-1,1]);%定坐标范围结果输入序列起点:n0=-6输入序列长度:N=15输入位移:n1=31011例7.3序列的移位和周期延拓已知，利用MATLAB生成并图示表示x(n)以8为周期的延拓）和，其中N=24。解：周期延拓方法1：利用矩阵乘法和冒号运算x=[1234];y=x'*ones(1,3);y1=(y(:))'P176方法2：采用求余函数mody=x(mod(n,M)+1)可实现对x(n)以M为周期的周期延拓。加1是因为MATLAB向量下标只能从1开始。)(8.0)(8nRnxn8N8(),(),(())(),((())xnxnmxnRnxn)())((N8nRmnx12%第七章例7.3程序q703序列x(n)的移位与周期延拓%N:观察窗口长度%M:序列x(n)长度%m:移位样点数clear;closeallN=24;M=8;m=input('输入移位值：m=');if(m1|m=N-M+1)%检验输入参数m是否合理error('输入数据不在规定范围内！');endn=0:N-1;x1=(0.8).^n;x2=[(n=0)&(nM)];xn=x1.*x2;%产生x(n)xm=zeros(1,N);%设定xm的长度fork=m+1:m+Mxm(k)=xn(k-m);end%xn移位m位，变为xmxc=xn(mod(n,8)+1);%产生x(n)的周期延拓xcm=xn(mod(n-m,8)+1);%产生x(n)周期延拓后的移位13subplot(2,2,1);stem(n,xn,'.');title('x(n)');axis([0,length(n),0,1]);ylabel('x(n)');xlabel('n');subplot(2,2,2);stem(n,xm,'.');title('x(n-m)');axis([0,length(n),0,1]);ylabel('x(n-m)');xlabel('n');subplot(2,2,3);stem(n,xc,'.');title('x(n)的周期延拓');axis([0,length(n),0,1]);ylabel('x(mod(n,8))');xlabel('n');subplot(2,2,4);stem(n,xcm,'.');title('x(n)的循环移位');axis([0,length(n),0,1]);ylabel('x(mod(n-m,8))');xlabel('n');结果：输入移位值：m=31415如果，xcm=xn(mod(n-m,8)+1);改为：xcm=xm(mod(n,8)+1);%产生x(n)移位后的周期延拓16例7.4离散系统对信号的响应本题给定6阶低通数字滤波器的系统函数，求它在下列输入序列x(n)下的输出序列y(n)。解：本题的计算原理见例6.14，在这里用工具箱函数filter来解。如果已知系统函数H(z)=B(z)/A(z)，则filter函数可求出系统对输入信号x(n)的响应y(n)。y=filter(B,A,x)P177例题7.4由差分方程可得到H(z)的分子和分母多项式系数向量A和B，再给出输入向量x即可。17%第七章例7.4程序q704求时域离散系统对常用序列的响应%N:输入数据长度%B:差分方程系数%A:差分方程系数clear;closeall;N=64;n=0:N-1;m=10;%下面设定系统参数A,BB=0.0003738*conv([1,1],conv([1,1],conv([1,1],conv([1,1],conv([1,1],[1,1])))))A=conv([1,-1.2686,0.7051],conv([1,-1.0106,0.3583],[1,-0.9044,0.2155]))x1=[n==0];%产生输入信号x1(n)y1=filter(B,A,x1);%对x1(n)的响应x2=[(n-m)==0];%产生输入信号x2(n)y2=filter(B,A,x2);%对x2(n)的响应x3=[n=0];%产生输入信号x3(n)y3=filter(B,A,x3);%对x3(n)的响应x4=[(n=0)&(n32)];%产生输入信号x4(n)y4=filter(B,A,x4);%对x4(n)的响应x5=exp(j*pi*n/8).*x4;%产生输入信号x5(n)y5=filter(B,A,x5);%对x5(n)的响应18subplot(3,2,1);stem(n,y1,'.');line([0,N],[0,0])axis([0,N,min(y1),max(y1)]);ylabel('y1(n)')subplot(3,2,2);stem(n,y2,'.');line([0,N],[0,0])axis([0,N,min(y2),max(y2)]);ylabel('y2(n)')subplot(3,2,3);stem(n,y3,'.');line([0,N],[0,0])axis([0,N,min(y3),max(y3)]);ylabel('y3(n)')subplot(3,2,4);stem(n,y4,'.');line([0,N],[0,0])axis([0,N,min(y4),max(y4)]);ylabel('y4(n)')subplot(3,2,5);stem(n,real(y5),'.');line([0,N],[0,0])axis([0,N,-1,1]);ylabel('Re[y5(n)]')subplot(3,2,6);stem(n,imag(y5),'.');line([0,N],[0,0])axis([0,N,-1,1]);ylabel('Im[y5(n)]')19运行结果：B=0.00040.00220.00560.00750.00560.00220.0004A=1.0000-3.18364.6223-3.77951.8136-0.48000.054420例7.5系统线性性质验证设系统差分方程为y(n)=x(n)+0.8y(n-1)，要求用程序验证系统的线性性质。解：产生两种输入序列，分别乘以常数后1.分别激励系统，再求输出之和；2.先相加，再激励系统求输出；对两个结果进行比较，方法是求它们之差，按误差的绝对值是否极小进行判断。21%Q705时域离散系统线性性质验证%N:输入信号长度%B:差分方程系数%A:差分方程系数clear;closeall;N=64;n=0:N-1;m=20;B=1;A=