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.*(n0);stem(t);%通过改变n.*(ni)或(n-i).*(ni)来实现不同的单位斜坡函数时移；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)%Implementsy(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)%Implementsy(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)%Implementsy(m+k)=x(n)n=m+n0;y=x;4.序列折叠：function[y,n]=sig_fold(x,m)%Implementsy(n)=x(-n)y=fliplr(x);n=-fliplr(n)5.序列奇偶性：function[xeven,xodd,m]=sig_even_odd(x,n)%Realsignaldecompositionintoevenandoddpartsif(imag(x)~=0)error('xisnotarealsequence');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)ifany(imag(x)~=0)error('Thesequenceisnotreal.')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);0510152025020400510152025020400510152025-200206.序列圆周移位：functiony=cirshift(x,m,N)%x:inputsequence%m:shiftnumber%N:showlengthiflength(x)Nerror('Nmustbegreaterthenlength(x)');endx=[xzeros(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)functionm=sigmod(n,N)m=rem(n,N);m=m+N;m=rem(m,N);例：由下例可以看出，循环移位不会改变频谱的幅值，只会改变其相位。n=0:10;x=10*0.8.^n;X=DFT(x,11);subplot(321);stem(x);subplot(323);stem(abs(X));subplot(325);stem(angle(X));y=cirshift(x,5,11);Y=DFT(y,11);subplot(322);stem(y);subplot(324);stem(abs(Y));subplot(326);stem(angle(Y));05101505100510150204060051015-10105101505100510150204060051015-5057.信号的循环/圆周卷积：functiony=circonvt(x1,x2,N)%计算循环卷积iflength(x1)Nerror('Length(x1)isnotgreaterthenN');endiflength(x2)Nerror('Length(x2)isnotgreaterthenN');endx1=[x1,zeros(1,N-length(x1))];x2=[x2,zeros(1,N-length(x2))];m=0:N-1;x2=x2(mod(-m,N)+1);H=zeros(N,N);forn=1:NH(n,:)=cirshift(x2,n-1,N);endy=x1*H';运行：x1=[122];x2=[1234];y=circonvt(x1,x2,5)y=9491414y=circonvt(x1,x2,4)y=1512914y=circonvt(x1,x2,7)y=149141480y=circonvt(x1,x2,10)y=149141480000以上方法是信号序列的时域循环卷积方法，当然也可以使用频域的卷积方法：先将两个序列进行DFT变换，再把DFT变换的结果进行点乘，然后对相乘的结果进行IDFT变换，即得到卷积结果移位：()()ynxnm；MATLAB实现：y=[zeros(1,m),x]折叠：y(n)=x(-n)；MATLAB实现：y=fliplr(x)；采样和：()()ynxn，MATLAB实现：y=sum(x(n1:n2))采样积：()()ynxn，MATLAB实现：y=prod(x(n1:n2))信号能量：2|()|xnExn,MATLAB实现：Ex=sum(abs(x).^2)