MATLAB在信号处理中的应用

zhouyang664

刚看了网友一个关于信号处理的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');

end

m=-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.')

end

N=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);

6.
序列圆周移位：

function y=cirshift(x,m,N)

%x:input sequence

%m:shift number

%N:show length

if length(x)>N

error('N must be greater then length(x)');

end

x=[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);

例：由下例可以看出，循环移位不会改变频谱的幅值，只会改变其相位。

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));

7.
信号的循环/圆周卷积：

function y=circonvt(x1,x2,N)

%计算循环卷积

if length(x1)>N

error('Length(x1) is not greater then N');

end

if length(x2)>N

error('Length(x2) is not greater then N');

end

x1=[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);

for n=1:N

H(n,:)=cirshift(x2,n-1,N);

end

y=x1*H';

运行：

>> x1=[1 2 2];x2=[1 2 3 4];

>> y=circonvt(x1,x2,5)

y =

9
4
9
14
14

>> y=circonvt(x1,x2,4)

y =

15
12
9
14

>> y=circonvt(x1,x2,7)

y =

1
4
9
14
14
8
0

>> y=circonvt(x1,x2,10)

y =

1
4
9
14
14
8
0
0

0
0

以上方法是信号序列的时域循环卷积方法，当然也可以使用频域的卷积方法：先将两个序列进行DFT变换，再把DFT变换的结果进行点乘，然后对相乘的结果进行IDFT变换，即得到卷积结果

移位： ；MATLAB实现：y=[zeros(1,m),x]

折叠：y(n)=x(-n)；MATLAB实现：y=fliplr(x)；

采样和： ，MATLAB实现：y=sum(x(n1:n2))

采样积： ，MATLAB实现：y=prod(x(n1:n2))

信号能量： ,MATLAB实现：Ex=sum(abs(x).^2)

因为有些公式是在MathType中编辑的，无法显示，这里附上源文件，参照着学习一下：

MATLAB在信号处理中的应用.doc (65.5 KB, 下载次数: 38)

2010-12-17 20:23 上传

点击文件名下载附件
MATLAB在信号处理中的应用

tenglang

thanks for sharing!
学习了

huzhi1988

谢谢咯。。。

zhq8846

很好的东西

ChaChing

顺著帖子上浮, 再次学了下!
忘记在那也看过sig_add/sig_mult这两个函数了, 细看下个人感觉好像可改进
在找序列的起始/终止下标, 使用min(n1), max(n1)...等, 好像根本多餘! min(n1)不就是n1(1), max(n1)不就是n1(end)! 难道时间序列可以不是单调递增吗?
所以练习下, 依个人习惯并无特别优化

1. function [f,n]=sigadd(f1,n1,f2,n2)
   
2. n=min(n1(1),n2(1)):max(n1(end),n2(end)); f=zeros(1,length(n));
   
3. ppp=n1(1)<=n & n<=n1(end); f(ppp)=f1;
   
4. ppp=n2(1)<=n & n<=n2(end); f(ppp)=f(ppp)+f2;

复制代码

动动熊

yalsrs

学习了。