MATLAB在信号处理中的应用

MATLAB在信号处理中的应用1、MATLAB在采样与波形发生的应用2、MATLAB在数字滤波器中的应用3、Simulink信号处理仿真的应用介绍1、MATLAB在采样与波形发生的应用数字信号处理的对象，是在采样时钟的控制之下，通过A/D转换器一定的采样频率对模拟信号进行采样得到的。根据采样定理，采样率必须大于模拟信号的最高采样率（Nyquist）的2倍。但是在许多情况下，要求信号以不同的频率采样，改变采样时钟虽然可行，但是并不可取。这时需要对采样数据进行处理：或者用抽取的方法降低其采样率（下采样）；或者用内插的方法提高其采样率（上采样）；或者两者兼有之（重采样）。在程序设计阶段，为了对程序进行调试或验证算法的正确性，需要一些特性已知的信号（简单的如正弦、方波、三角波等），所以可以由Matlab提供一些波形产生函数。1、MATLAB在采样与波形发生的应用常用的信号发生函数有：（1）upfirdn（）改变信号的采样率（用FIR滤波器）（2）decimat（）下采样（先经低通滤波）（3）interp（）上采样（使用了的通滤波器）（4）resample（）改变信号的采样率（5）chirp（）线性调频信号发生器（6）diric（）周期sinc函数（Dirichlet函数）（7）gauspuls（）高斯函数调幅的正弦波发生器（8）pulstran（）通过对连续函数或脉冲原型进行采样而得到的脉冲序列的发生器（9）rectpuls（）非周期矩形波发生器（10）sawtooth（）锯齿波和三角波发生器（11）sinc（）sinc函数发生器（12）square（）方波发生器（13）tripuls（）非周期三角脉冲发生器下面介绍一些常见的信号函数：1、MATLAB在采样与波形发生的应用1.1连续信号函数1、指数信号形式：y=A*exp（a*t）2、正/余弦信号形式：y=A*sin（w0*t+f）

或：y=A*cos（w0*t+f）3、抽样函数用sinc函数表示，其定义为：sinc(t)=sin(πt)/(πt)MATLAB语言调用形式:y=sinc(t)1、MATLAB在采样与波形发生的应用举例：t=-3*pi:pi/100:3*pi;y=sinc(t/pi);plot(t,y),gridon1、MATLAB在采样与波形发生的应用4、矩形脉冲信号MATLAB语言调用形式:y=rectpuls(t,width)产生一个幅值为1,宽度为width、相对于t=0点左右对称的矩形波信号t=0:0.001:4;T=1;y=rectpuls(t-2*T,2*T);plot(t,y),gridonaxis([0,4,-0.5,1.5])%以t-2*T=0,即t=2*T为对称中心1、MATLAB在采样与波形发生的应用周期性矩形脉冲信号MATLAB语言调用形式:y=square(t,duty)产生一个周期为2π,幅值为±1的周期性方波,duty表示占空比(dutycycle)t=-0.1:0.001:0.1;y=2*square(2*pi*30*t,75);plot(t,y)gridonaxis([-0.1,0.1,-2.5,2.5])%w=2*pi*30%占空比75%1、MATLAB在采样与波形发生的应用5、三角脉冲信号MATLAB语言调用形式:y=tripuls(t,width,skew)产生一个最大幅值为1,宽度为width、斜度为skew的三角波信号.t=-3:0.001:3;y=tripuls(t,4,0.5);plot(t,y)gridonaxis([-3,3,-1.5,1.5])1、MATLAB在采样与波形发生的应用周期性三角脉冲信号MATLAB语言调用形式:y=sawtooth(t,width)产生一个周期为2π,最大幅值为1,最小幅值为-1周期性三角波信号.t=-5*pi:0.0001:5*pi;y=sawtooth(t,0.5);plot(t,y)gridonaxis([-16,16,-1.5,1.5])1、MATLAB在采样与波形发生的应用6、一般周期性脉冲信号MATLAB语言调用形式:y=pulstran(t,d,’func’)依据名为func的连续函数并以之为一个周期,从而产生一串周期性的连续函数。d用于指定周期性的偏移量，这个func函数会被计算length(d)次,最后值相当于:y=func(t-d(1))+func(t-d(2))+...一般的调用形式为:y=pulstran(t,d,’func’,p1,p2)其中，p1,p2为func函数需要的额外输入参数值,最后值相当于:y=func(t-d(1),p1,p2,...)+func(t-d(2),p1,p2,...)+程序举例：t=0:0.0001:1;d=0:1/3:1;y=pulstran(t,d,'rectpuls',0.1);figure(1)plot(t,y)gridonaxis([0,1,-0.1,1.5])%p1=0.1，为带宽yy=pulstran(t,d,'tripuls',0.1,-1);figure(2)plot(t,yy)gridonaxis([0,1,-0.1,1.5])%0.1,1位宽度和斜度1、MATLAB在采样与波形发生的应用举例：产生一频率为10kHz的周期高斯脉冲信号，其带宽为50％。脉冲重复的频率为1kHz，采样率为50kHz，脉冲序列的长度为10ms。重复时幅度每次衰减为原来的0.8倍。t=0:1/50e3:10e-3;d=[0:0.001:0.01;0.8.^(0:10)]';y=pulstran(t,d,'gauspuls',1000,0.5);plot(t,y);xlabel('时间/s');ylabel('幅值')1、MATLAB在采样与波形发生的应用1.2离散信号函数1、指数序列Y=a.^k2、正弦序列离散的指数序列与正弦序列用MATLAB的表示与连续信号类似,只不过用stem函数而不是用plot函数来画出序列的波形.3、单位脉冲序列1、MATLAB在采样与波形发生的应用k=-50:50;y=[zeros(1,50),1,zeros(1,50)];stem(k,y)1、MATLAB在采样与波形发生的应用4、单位阶跃序列k=-50:50;y=[zeros(1,50),ones(1,51)];stem(k,y),axis([-60,60,-0.5,1.5])5、变换函数

1）一维快速离散Fourier变换：语法格式：y=fft(x)。y是计算信号x的快速离散傅里叶变换。当x为矩阵时，计算x中每一列信号的离散傅里叶变换。当x的长度为2的幂时，用基2算法；否则，采用较慢的分裂基算法。

y=fft(x,n)。计算n点的FFT。当x的长度大于n时，截断x;否则补零。2）一维快速离散Fourier逆变换：

语法格式：y=ifft(x)。y是计算信号x的快速离散傅里叶变换逆变换。

y=ifft(x,n)。计算n点的快速离散傅里叶变换的逆变换。3）离散余弦变换（DCT）

语法格式：y=dct(x)。计算信号x的离散余弦变换。

y=dct(x,n)。计算n点的离散余弦变换。当x的长度大于n时，截断x；否则补零。离散余弦逆变换可由函数idct实现。2、MATLAB在数字滤波器中的应用滤波器是指用来对输入信号进行滤波的硬件或软件。如果滤波器的输入、输出都是离散时间信号，则该滤波器的冲击响应也必然是离散的，这样的滤波器定义为数字滤波器。数字滤波器用硬件实现的基本部件包括延迟器、乘法器和加法器；如果用软件来实现时，它即是一段线性卷积程序。软件实现的优点是系统函数具有可变性，仅依赖于算法结构，而且易于获得较理想的滤波功能。Matlab的信号处理工具箱的两个基本组成就是滤波器的设计与实现以及谱分析。工具箱提供了丰富而简便的设计、实现FIR和IIR的方法，使原来繁琐的程序设计简化成函数调用，特别是滤波器的表达方式和形式之间的相互转换显得十分简便，为滤波器的设计和实现开辟了一片广阔的天地。2、MATLAB在数字滤波器中的应用2.1、模拟低通滤波器MATLAB中给出了几种不同类型的模拟滤波器原型。低通模拟Butterworth滤波器原型——buttap

[Z,P,K]=buttap(n)；

其中，Z，P，K分别是n阶Butterworth滤波器的零点、极点和增益。低通模拟ChebyshevI型滤波器原型——cheb1ap

[Z,P,K]=cheb1ap(n,rp)其中，Z，P，K分别是n阶ChebyshevI型滤波器的零点、极点和增益。滤波器在通带内的最大衰减为rp。ChebyshevI型滤波器的主要特点是在阻带内达到最大平滑。低通模拟ChebyshevⅡ型滤波器原型——cheb2ap

[Z,P,K]=cheb2ap(n,rs)2、MATLAB在数字滤波器中的应用其中，Z，P，K分别是n阶ChebyshevⅡ型滤波器的零点、极点和增益。滤波器在阻带内的最小衰减为rs。ChebyshevⅡ型滤波器的主要特点是在通带内达到最大平滑。低通模拟椭圆滤波器原型——ellipap

[Z,P,K]=ellipap(n,rp,rs)其中，

Z，P，K分别是n阶椭圆滤波器的零点、极点和增益。滤波器在通带内的最大衰减为rp，在阻带内的最小衰减为rs。低通模拟Bessel滤波器——besselap

[Z,P,K]=besselap(n)其中，

Z，P，K分别是n阶低通模拟Bessel滤波器的零点、极点和增益。2、MATLAB在数字滤波器中的应用分别设计20阶ChebyshevI型低通模拟滤波器，通带内的最大衰减为0.3dB；20阶低通模拟ChebyshevⅡ型滤波器阻带内的最小衰减为45dB，并给出其频率特性图。[z1,p1,k1]=cheb1ap(20,0.3);[num1,den1]=zp2tf(z1,p1,k1);%将零极点形式转换为系统函数形式[z2,p2,k2]=cheb2ap(20,45);[num2,den2]=zp2tf(z2,p2,k2);figure(1)freqs(num1,den1)%求频率响应

figure(2)freqs(num2,den2)2、MATLAB在数字滤波器中的应用2、MATLAB在数字滤波器中的应用对于模拟高通、带通、带阻滤波器，其设计方法为先将要设计的滤波器的技术指标通过某种频率转换关系转换成低通滤波器的技术指标，并依据这些指标设计出低通滤波器的转移函数，然后再依据频率转换关系变成所设计的滤波器的转移函数。

MATLAB的信号处理工具箱提供了从低通滤波器向低通、高通、带通、带阻滤波器转换的函数。（1）低通向低通的转变——lp2lp

[numt,dent]=lp2lp(num,den,Wo)

[At,Bt,Ct,Dt]=lp2lp(A,B,C,D,Wo)

这两条语句分别将表示成传递函数形式与状态方程形式的具有任意截止频率的模拟低通滤波器原型转换成截止频率为Wo的低通滤波器。2、MATLAB在数字滤波器中的应用（2）低通向高通的转变——lp2hp

[numt,dent]=lp2hp(num,den,Wo)

[At,Bt,Ct,Dt]=lp2hp(A,B,C,D,Wo)

这两条语句分别将表示成传递函数形式与状态方程形式的具有任意截止频率的模拟低通滤波器原型转换成截止频率为Wo的高通滤波器。3）低通向带通的转变——lp2bp

[numt,dent]=lp2bp(num,den,Wo,Bw)

[At,Bt,Ct,Dt]=lp2bp(A,B,C,D,Wo,Bw)这两条命令分别将表示成传递函数形式与状态方程形式的具有任意截止频率的模拟低通滤波器原型转换成中心频率为Wo、带宽为Bw的带通滤波器。（4）低通向高通的转变——lp2bs

[numt,dent]=lp2bs(num,den,Wo,Bw)

[At,Bt,Ct,Dt]=lp2bs(A,B,C,D,Wo,Bw)

这两条命令分别将表示成传递函数形式与状态方程形式的具有任意截止频率的模拟低通滤波器原型转换成中心频率为Wo、带宽为Bw的带阻滤波器2、MATLAB在数字滤波器中的应用采用两种方法设计10阶模拟椭圆低通滤波器，通带的最大衰减为3dB，阻带内的最大衰减为40dB，截止频率为8πrad。[z,p,k]=ellipap(10,3,40);

%根据要求设计零极点形式低通模拟椭圆滤波器[A,B,C,D]=zp2ss(z,p,k);

%零极点形式向状态方程形式转换[at,bt,ct,dt]=lp2lp(A,B,C,D,8*pi);

%给定截止频率的低通滤波器[num2,den2]=ss2tf(at,bt,ct,dt)figurefreqs(num2,den2)%求频率响应

%%第二种方式[z2,p2,k2]=ellipap(10,3,40);[num,den]=zp2tf(z,p,k);%零极点形式向传递函数形式转换[num1,den1]=lp2lp(num,den,8*pi);%给定截止频率的低通滤波器figurefreqs(num1,den1)2、MATLAB在数字滤波器中的应用2、MATLAB在数字滤波器中的应用为了设计数字滤波器，需要将模拟滤波器离散化，两个离散化基本方法为冲激响应不变法和双线性变换法，MATLAB中也有两种相应的函数。（1）脉冲响应不变法脉冲响应不变法的基本原理：从滤波器的脉冲响应g(t)出发，对其以周期T进行采样得到离散序列g(nT)作为数字滤波器的脉冲响应。相应的MATLAB函数为impinvar，具体格式为

[bz,az]=impinvar(b,a,fs)

将具有[b,a]模拟滤波器传递函数模型转换成采样频率为fs的数字滤波器的传递函数[bz,az]。如果在函数中没有指定采样频率fs，默认为1Hz。2、MATLAB在数字滤波器中的应用（2）双线性变换法为了克服脉冲响应不变法产生的频率混叠现象，模拟向数字的转变通常采用双线性变换法。相应的MATLAB函数为bilinear，具体格式为（1）[zd,pd,kd]=bilinear(z,p,k,fs)（2）[numd,dend]=bilinear(num,den,fs)（3）[Ad,Bd,Cd,Dd]=bilinear(A,B,C,D,fs)

以上三个命令分别将零极点形式、传递函数形式和状态方程形式的模拟滤波器转换成相应形式的数字滤波器，其中fs是采样频率。2、MATLAB在数字滤波器中的应用2.2、无限长单位脉冲响应滤波器–IIRMATLAB设计IIR数字滤波器可分为以下几步来实现：（1）按一定规则将数字滤波器的技术指标转换成模拟低通滤波器的技术指标；（2）根据转换后的技术指标确定滤波器的最小阶数和截止频率；（3）利用最小阶数产生模拟低通滤波器原型；（4）利用截止频率把低通滤波器原型转换成模拟低通、高通、带通或带阻滤波器；（5）利用冲激响应不变法或双线性变换法把模拟滤波器转换成数字滤波器。2、MATLAB在数字滤波器中的应用函数名功能格式buttord计算Butterworth滤波器的阶数和截止频率[N,Wc]=buttord(wp,ws,Rp,Rs,’s’)cheb1ord计算ChebyshevI型滤波器的阶数[N,Wc]=cheb1ord(wp,ws,Rp,Rs,’s’)cheb2ord计算ChebyshevⅡ型滤波器的阶数[N,Wc]=cheb1ord(wp,ws,Rp,Rs,’s’)ellipord计算椭圆滤波器的最小阶数[N,Wc]=ellipord(wp,ws,Rp,Rs,’s’)其中，N为滤波器阶数，Wc为截止频率，wp、ws是通带和阻带的截止频率，Wc、wp、ws均为被采样频率归一化后的频率，范围是[0,1]。用到相关函数见表格：2、MATLAB在数字滤波器中的应用MATLAB还提供了几个直接设计IIR数字滤波器的函数，这些函数也要与数字滤波器的阶数选择配合使用。

如设计Butterworth数字滤波器的butter函数：（1）[b,a]=butter(n,Wn)：设计n阶截止频率为Wn的Butterworth低通数字滤波器的传递函数模型系数[b,a]，系数长度为n+1。Wn为归一化截止频率，其最大值为采样频率的一半。当Wn为两元素向量Wn=[w1,w2]时，函数返回2N阶的带通数字滤波器，通带为w1<w<w2。（2）[b,a]=butter(n,Wn,’high’)：设计高通滤波器系数b，a。（3）[b,a]=butter(n,Wn,’stop’)：设计带阻滤波器系数b，a，频率Wn=[w1,w2]，阻带为w1<w<w2。（4）[A,B,C,D]=butter(n,Wn)和[z,p,k]=butter(n,Wn)：分别返回所设计数字滤波器的状态模型系数[A,B,C,D]和零极点增益[z,p,k]。2、MATLAB在数字滤波器中的应用一个数字滤波器的抽样频率为fs=2kHz，试设计一个为此系统使用的带通滤波器，希望采用椭圆型滤波器。要求（1）通带范围为300Hz到400Hz，在通带边缘频率处的衰减不大于0.5dB，（2）在200Hz以下和500Hz以上衰减不小于15dB。fs=2000;ws=[300,400]/(fs/2);wp=[200,500]/(fs/2);rp=0.5;rs=15;

%滤波器参数转换[N,wc]=ellipord(wp,ws,rp,rs);%阶数与截止频率[num,den]=ellip(N,rp,rs,wc);%直接设计椭圆滤波器[H,W]=freqz(num,den);F=W*fs/2/pi;plot(F,20*log10(abs(H)));gridxlabel('频率（Hz）');ylabel('幅值（dB）');axis([01000-500])2、MATLAB在数字滤波器中的应用2、MATLAB在数字滤波器中的应用2.3、有限长单位脉冲响应滤波器—FIR有限长单位脉冲响应滤波器FIR的系统函数为：

FIR系统只有零点，因此FIR系统不像IIR系统那样易取得较好的通带和阻带衰减特性。但FIR系统具有自己突出的优点：稳定且具有线性相位。2、MATLAB在数字滤波器中的应用有限长单位脉冲响应FIR滤波器的设计方法建立在对理想滤波器特性作某种近似的基础上，主要包括窗函数法和频率采样法。（1）.窗函数法窗函数法设计的原理是，用一个截短后的单位抽样响应h(n)来逼近理想的非因果的单位抽样响应hd(n)。即

上式为窗函数。MATLAB信号处理工具箱提供了Boxcar（矩形窗）、Bartlet（巴特利特）、Hanning（汉宁）等窗函数2、MATLAB在数字滤波器中的应用矩形窗：

w=boxcar(M)：返回M点矩形窗序列。其频率响应的逼近程度好坏完全取决于窗函数的频率特性。除了提供窗函数命令外，MATLAB还提供了专用窗函数法设计FIR滤波器的命令fir1。

b=fir1(n,Wn,'ftype')：用来设计具有线性相位的n阶FIR数字滤波器，返回的向量b为滤波器的系数。截止频率Wn必须在0到1之间，数值1对应于fs/2。滤波器的归一化增益在Wn处为-6dB。当Wn=[W1W2]时，fir1返回一个N阶的带通滤波器，通带为[W1W2]。2、MATLAB在数字滤波器中的应用用矩形窗设计线性相位FIR低通滤波器，通带截止频率为0.2，滤波器单位脉冲响应的长度为21。M=21;wc=0.2*pi;r=(M-1)/2;nr=-r:r;hdn=sin(wc*nr)/pi./(nr+eps);ifrem(M,2)~=0,hdn(r+1)=wc/pi;endw=boxcar(M);%矩形窗hn=hdn.*w‘;%时域脉冲响应hw=fft(hn,512);w1=2*[0:511]/512;%频谱特性subplot(2,1,1),stem(1:M,hn,‘o’);%针状图xlabel('n'),ylabel('h(n)');subplot(2,1,2),plot(w1,20*log10(abs(hw)));xlabel('w/pi');ylabel('幅值（dB）');axis([02-9010]);2、MATLAB在数字滤波器中的应用2、MATLAB在数字滤波器中的应用（2）频率采样法

频域采样法的基本原理是对理想的滤波器的频率响应在频域进行采样，以此来确定FIR滤波器的H(k)。在MATLAB信号处理箱中，为频率采样法设计FIR滤波器提供了专用函数命令fir2。该函数所得FIR数字滤波器系数为实数，具有线性相位，且满足对称性。2、MATLAB在数字滤波器中的应用fir2命令的格式如下：1）b=fir2(n,f,m)：返回n阶FIR数字滤波器的系数b。向量f和m分别指定滤波器采样点的频率和幅值，f为归一化频率。所期望的滤波器的频率响应可用plot(f,m)绘制。2）b=fir2(n,f,m,window)：用指定的窗函数来设计FIR数字滤波器，窗函数包括Boxar、Hann、Bartlett、Blackman、Kaiser及Chebwin等。默认情况下，使用hamming窗2、MATLAB在数字滤波器中的应用例：使用fir2命令设计上例中的滤波器。f=0:0.01:1m(1:21)=1;m(22:41)=0;m(42:61)=1;m(62:81)=0;m(82:101)=1;%设置多通带plot(f,m,'k:');holdonN=30;b=fir2(N,f,m);%令第一个频带为通带[h,f1]=freqz(b);f1=f1./pi;plot(f1,abs(h));legend('理想滤波器','设计滤波器');xlabel('归一化频率'),ylabel('幅值');2、MATLAB在数字滤波器中的应用2、MATLAB在数字滤波器中的应用下面介绍其它几种信号检测方法在MATLAB中的应用1、相关检测法1.1、自相关检测基本原理：利用信号周期性和噪声随机性的特点，通过自相关或互相关运算，打到除去噪声的目的。自相关函数：[a,b]=autocorr(x,'option');互相关函数：[c,lags]=xcorr(x,lag,'option')；其中，lags为延迟矢量，其范围为[-maxlags:maxlags]'option'='unbiased'or'biased';2、MATLAB在数字滤波器中的应用2、互相关检测基本原理：若已知发送信号的重复周期或频率，就可以在接收端发出一个与发送信号周期相同的“干净的”重复周期信号，称为本地信号。将它与混有噪声的信号进行互相关，即可除去噪声的影响，提高电路的抗干扰能力。应用：锁定放大器，可以把待测信号中与参考信号同步的信号放大并检测出来。函数：[a,b]=xcorr(x,y,'unbiased');上式与自相关检测相似，多添加的y可表示x与y的互相关检测3、小波（包）分析“小波”就是小的波形，所谓“小”是指它具有衰减性；而称之为“波”则是指它的波动性，其振幅正负相间的震荡形式。小波变换是时间（空间）频率的局部化分析，它通过伸缩平移运算对信号(函数)逐步进行多尺度细化，最终达到高频处时间细分，低频处频率细分，能自动适应时频信号分析的要求，从而可聚焦到信号的任意细节。信号处理的目的就是：准确的分析、诊断、编码压缩和量化、快速传递或存储、精确地重构（或恢复）。从数学地角度来看，信号与图像处理可以统一看作是信号处