matlab基础4

1.fft ()：计算函数的快速傅里叶变换。

FFT是DFT的快速算法。

调用格式：Y = fft(x)：直接计算x(n)的DFT，计算长度为序列x(n)的长度Y = fft(x,n)：如果x(n)的长度小于n，则在其后补0，使之成为长度为n的序列，即计算n点DFT，结果产生n个数据。

注：1)fft函数是由机器语言而不是Matlab指令写成的，是Matlab的一个内建函数。

2)同样逆变换也是由内建函数ifft实现的。

3)Fftshift()调用格式：Y= fftshift(X)，用来重新排列X= fft(x)的输出。

即Y序列是X 序列的移位序列，即将X的零频分量移到频谱的中心。

4)Abs()：（求绝对值）计算信号的幅频特性。

5)Angle():计算信号的相频特性。

6)Floor()：向下取整函数k=floor(-(N-1)/2:(N-1)/2);表示频序向量值。

2. 在各种傅利叶变换中，能够用计算机以高效计算的，只有离散傅利叶变换DFT ，所以要计算任意信号的频谱，必须弄清各种变换与DFT 的关系。

⑴非周期连续信号的频谱计算。

设： 为绝对可积的连续信号，则其数字频谱 可以用FFT 求得。

即用下式可求得连续信号的频谱。

第一步：先讨论信号时域的长度L 。

例如; ，此信号的峰会幅度为1，而在t>50时，其幅度小于0.0067，因些如果选择L=50，记录时间长度将覆盖信号的主要部分。

()()[]n x DFT T j X k ⋅≈ΩΩ=Ω0()t x a ()()Ω=j X X ω()t a e t x 1.0-=第二步：因为这个信号不是有限带宽，则时域采样所造成的频率混叠可以忽略。

第三步：频率分辨率。

原连续信号的频谱是定义在所有的频率点上，但使用FFT 来计算时，只计算出了它的N 个频率值。

所以必须从这N 个值中进行插值运算求原连续的频谱。

所以频率的分辨率越小，插值的效果越好。

而这个分辨率与长度L 成反比，因为所以我们愿意选择较大的L 来改善频率分辨率。

第四步：考虑截断效应。

由于对无限长信号进行截断，这样就会把窗函数频谱引入信号频谱，引起吉布斯效应（波动），同时也产生频谱混叠。

综上所述：先初步选择时间记录长度L ，使得在此之间包括大部分非零的值。

然后用逐渐减小T 和加大N 的步骤来选择周期T 。

方法是：选择一个较小的N 值，算出T=L/N ，然后用这个N 值求出一个频谱。

再选择一个N 值，重复这一计算，比较结果。

因为它们是在同一个L 下产生的，所以它们具有相同的截断效应，这样两个频谱间的误差就是由频率混叠引起的。

如果两个频谱非常接近，则频率混叠可以忽略。

如果精度不够，可以另选择较小的T ，再进行重复计算。

()⎪⎭⎫ ⎝⎛=⎪⎭⎫ ⎝⎛=⎪⎭⎫ ⎝⎛=L k X NT k X D N k X kD X πππ222选定T 之后，再进一步选择L 。

例如可将其增加一倍，再与前面的结果比较。

注意：可用乃奎斯特频率处的幅频特性来评价混叠的严重程度。

⑵连续周期信号的频谱计算。

例如 ，它的理论频谱是 包含了权重为pi 的位于 上的两个冲击。

理论上，只要采样周期小于5π=0.6283s ，就不会发生频率混叠。

但是由于将它截断了，所以频谱不再是有限带宽，即必须采用更小的采样周期。

()t t x a 5cos =()()[]()[]55+Ω+-Ω=ΩδπδπX 5±=Ω3.结论：用DFT去计算一个非周期信号的傅利叶变换，则其频谱的正常幅度与用DFT计算相差一个加权T。

T=0.002; %采样间隔T=0.002F=1/T; %采样频率为F=1/TL=0.1; %记录长度L=0.01，体现所研究模拟信号的长度N=L/T;t=0:T:L;a=2;f=0:F/N:F; %f=k/NT,k为整数，k从0到Nx=exp(-10*t);y=T*abs(fft(x)); %有限长（长度为N)离散时间信号x1的dft再乘T来近似模拟信号的频谱，长度为Nsubplot(2,1,1),stem(t,x);title('指数信号');subplot(2,1,2),stem(f,y);title('指数信号频谱');4.用三种方法进行两序列的圆周卷积：时域手工方法调用Matlab函数function y = circonv(x,h,N) 用频域DFT相乘，再求反变换。

注：用函数：L=length(x)，计算序列x的长度。

5.滤波器函数介绍系统函数的负幂形式：()MMMMzazazbzbbzH----++++++=11111，也被称为传递函数。

系统函数的零极增益形式：()()()∏∏==---=NkkMllMNpzzzzbzH11IIR型滤波器的设计流程：1.设计模拟低通滤波器：如巴特沃斯低通滤波器，MATLAB提供函数N阶归一化模拟原型滤波器。

其调用格式为[z,p,k] = buttap(n)；只要输入阶数N，就可返回零点数组z和极点数组p以及增益k。

2.阶数N由函数3.进行频带变换，从归一化模拟低通变换到具有特定截止频率的模拟低通、高通、带通、带阻滤波器。

MATLAB提供lp2lp lp2hp lp2bp lb2bs这四个函数进行转换。

但是这个函数调用的参数是系统函数的负幂形式中的分子，分母系数。

所以需要将利用函数[z,p,k] = buttap(n)求出的零极增益形式转换成负幂形式。

MATLAB提供函数进行转换。

4.利用双线性变换（或冲击响应不变法）将模拟滤波器转换为数字滤波器。

MATLAB提供函数5.下面利用函数()ωj e H。

用6.用设计好的滤波器对已知序列进行滤波。

MATLAB提供函数Butterworth filter order and cutoff frequency[N,Wn] = buttord(Wp,Ws,Rp,Rs)[N,Wn] = buttord(Wp,Ws,Rp,Rs,'s')可根据给定的滤波器指标Wp,Ws,Rp,Rs，即通带截止频率，阻带截止频率，通带允许最大衰减，阻带允许最小衰减，求得滤波器的阶阶数N。

Wn为巴特沃斯滤波器固有的频率。

s说明所设计的是模拟滤波器。

Butterworth analog lowpass filter prototype这个函数来设计N阶归一化巴特沃斯模拟原型滤波器，其调用格式：为[z,p,k] = buttap(n)只要输入阶数N ，就可返回零点数组z 和极点数组p 以及增益k 。

Convert zero-pole-gain filter parameters to transfer function form 调用格式：[b,a] = zp2tf(z,p,k)它属于线性系统变换函数。

由传递函数变换为零极增益的函数：[z,p,k] = tf2zp(b,a)Change cutoff frequency for lowpass analog filter调用格式：[b,a]=lp2lp(Bap,Aap,Wn)MATLAB 信号处理工具箱提供了从归一化低通滤波器到低通、高通、带通和带阻滤波器的变换函数，只要已经设计好了一个归一化模拟低通滤波器，知道了它们的分子分母系数向量，就可以直接调用这些函数。

共有四个：lp2lp lp2hp lp2bp lb2bs 。

Wn 为截止频率。

[b,a]为按1-z 的升幂排列的多项式的分子和分母的系数向量，[Bap,Aap]也是如此。

()M M MM z a z a z b z b b z X ----++++++= 111101Bilinear transformation method for analog-to-digital filter conversion调用格式：[numd,dend]=bilinear(b,a,fs)MATLAB信号处理工具箱提供一个名为Bilinear的子程序，它用双线性变换法来实现把模拟滤波器变换为数字滤波器。

Fs为采样频率。

调用freqz函数的方式如下：[H,w] = freqz(b,a,N)，或N缺省。

输出变元H是系统的频率特性，它是一个N元的复数向量。

π输出变元w是数字频率向量，它把0到pi均分为N份，分辨率为N输入变元b和a分别为分子分母多项式的负幂系数向量，即多项式的首项应是常数项，以后按jwe-的升幂排列，由此形成多项式的系数向量。

输入变元N为所选的频率点数，它决定了频率分辨率的密度，N的默认值为512。

这样求出的频率特性是在正频率区间的特性。

调用格式：y = filter(b,a,X)求系统响应时，有一种方法是递推法，差分方程的递推过程在MATLAB 中用一个称为filter的子程序来实现，调用的最简单形式为：y = filter(b,a,X)。

其中[]bM b b b ,1,0=，[]aN a a a ,1,0=，是差分方程()()∑∑==-=-MI i Nk ki n x b k n y a的系数数组。

而X 则是输入序列。

X 通过一个系统的响应。

利用窗函数法设计FIR 型滤波器的流程：1.按设计要求确定理想的)(ωj d e H 计算出()n h d ，例如低通滤波器单位冲击响应：()()[]()απαω--=n n n h c d sin2. 根据过渡带宽，阻带最小衰减，确定窗函数的形状、长度。

例如：给定抽样频率为Ωst=3π×104(rad/s),通带截止频率为Ωp=3π×103(rad/s),阻带起始频率为Ωs=6π×103(rad/s),阻带衰减比小于50dB 。

计算：第一步：求理想低通滤波器的截止频率c ωπω3.02=Ω+Ω=Ω=Ω=Ω=stps st c st c c c f f f T π2stst f Ω=求出c ω，则可计算()()[]()απαω--=n n n h c d sin ，其中21-=N τ，是保证线性相位所必须的移位。

第二步：由于阻带衰减比小于50dB ，所以根据表可选择海明窗Hamming(N)或布莱克曼窗Blackman(N)。

若选择海明窗，则其过渡带宽为：NN A πω6.6==∆p s Ω-Ω=∆Ωst f T ∆Ω=∆Ω=∆ωππω2.06.6=∆Ω===∆∴stf N N A332.06.6==∴ππN1621=-=N τ 这样，我们确定先选择的窗函数为海明窗，其长度为３３3. 求所设计的FIR 滤波器的单位抽样响应 ，就完成了FIR 滤波器的设计。

在MATLAB 中，用语句：h=hd.*w1'; 可实现其中，hd 是理想的单位冲击响应，调用函数hd=ideallp(wc,N) 实现，其中wc 为理想低通的截止频率,N 为窗函数长度。