特别说明：完整45期数字信号处理教程，原创高性能示波器代码全开源地址：链接
; P- [/ _5 W+ D) q

第36章 FIR滤波器的Matlab设计

' P$ i8 X) S/ @9 M" e

    本章节讲解FIR滤波器的Matlab设计。主要是函数fir1和fir2的使用。

    36.1 窗函数

    36.2 fir1函数

    36.3 fir2函数

    36.4 总结

. B' v# v- A! C  y; U7 R* e  o36.1 窗函数

    在数字信号处理中不可避免地要用到数据截取的问题。例如，在应用DFT的时候，数据x(n)总是有限长的，在滤波器设计中遇到了对理想滤波器抽样响应h(n)的截取问题，在功率谱估计中也要遇到对自相关函数的截取问题。总之，我们在实际工作中所能处理的离散序列总是有限长，把一个长序列变换成有限长的序列不可避免的要用到窗函数。因此，窗函数本身的研究及其应用是信号处理中的一个基本问题。

    不同的窗函数对信号频谱的影响是不一样的，这主要是因为不同的窗函数，产生泄漏的大小不一样，频率分辨能力也不一样。信号的截断产生了能量泄漏，而用FFT算法计算频谱又产生了栅栏效应，从原理上讲这两种误差都是不能消除的，但是我们可以通过选择不同的窗函数对它们的影响进行抑制。(矩形窗主瓣窄，旁瓣大，频率识别精度最高，幅值识别精度最低；布莱克曼窗主瓣宽，旁瓣小，频率识别精度最低，但幅值识别精度最高)。

    对于窗函数的选择，应考虑被分析信号的性质与处理要求。如果仅要求精确读出主瓣频率，而不考虑幅值精度，则可选用主瓣宽度比较窄而便于分辨的矩形窗，例如测量物体的自振频率等；如果分析窄带信号，且有较强的干扰噪声，则应选用旁瓣幅度小的窗函数，如汉宁窗、三角窗等；对于随时间按指数衰减的函数，可采用指数窗来提高信噪比。

l 矩形窗：

    矩形窗属于时间变量的零次幂窗。矩形窗使用最多，习惯上不加窗就是使信号通过了矩形窗。这种窗的优点是主瓣比较集中，缺点是旁瓣较高，并有负旁瓣，导致变换中带进了高频干扰和泄漏，甚至出现负谱现象。

l 三角窗：

    三角窗亦称费杰（Fejer）窗，是幂窗的一次方形式。与矩形窗比较，主瓣宽约等于矩形窗的两倍，但旁瓣小，而且无负旁瓣。

l 汉宁窗：

    汉宁窗又称升余弦窗，汉宁窗可以看作是3个矩形时间窗的频谱之和，或者说是 3个 sinc（t）型函数之和，而括号中的两项相对于第一个谱窗向左、右各移动了 π/T，从而使旁瓣互相抵消，消去高频干扰和漏能。可以看出，汉宁窗主瓣加宽并降低，旁瓣则显著减小，从减小泄漏观点出发，汉宁窗优于矩形窗．但汉宁窗主瓣加宽，相当于分析带宽加宽，频率分辨力下降。

l 海明窗：

    海明窗也是余弦窗的一种，又称改进的升余弦窗。海明窗与汉宁窗都是余弦窗，只是加权系数不同。海明窗加权的系数能使旁瓣达到更小。分析表明，海明窗的第一旁瓣衰减为一42dB．海明窗的频谱也是由3个矩形时窗的频谱合成，但其旁瓣衰减速度为20dB/（10oct），这比汉宁窗衰减速度慢。海明窗与汉宁窗都是很有用的窗函数。

l 高斯窗：

    三角窗亦称费杰（Fejer）窗，是幂窗的一次方形式。与矩形窗比较，主瓣宽约等于矩形窗的两倍，但旁瓣小，而且无负旁瓣。

    还有很多其它的窗口这里就不做介绍了，需更详细的了解的话，可以看matlab中help文档中的如下部分（这里只截图了一部分窗口）。

    或者直接在命令窗口输入wintool可以打开窗口工具：

    打开后界面如下：

36.2 fir1函数
/ F4 F, M! w4 i2 |36.2.1 fir1函数介绍

    函数fir1用来设计标准频率响应的基于窗函数的FIR滤波器，可实现加窗线性相位FIR滤波器设计。

语法：

    b = fir1(n,Wn)

    b = fir1(n,Wn,'ftype')

    b = fir1(n,Wn,window)

    b = fir1(n,Wn,'ftype',window)

    b = fir1(...,'normalization')

其中，n：为了滤波器的阶数；Wn：为滤波器的截止频率；ftype：参数用来决定滤波器的类型，当ftype=high时，可设计高通滤波器，当ftype=stop时，可设计带阻滤波器。Window参数用来指导滤波器采用的窗函数类型。其默认值为汉明（Hamming）窗。

    使用fir1函数可设计标准的低通，高通，带通和带阻滤波器。滤波器的系数包含在返回值b中，可表示为：

        b(z) = b(1) + b(2)z-1 + …… +b(n+1)z-n

(1) 采用汉明窗设计低通FIR滤波器

使用b=fir1(n, Wn)可得到低通滤波器。其中， 0Wn1, Wn=1相当于0.5。其语法格式为

b=fir1(n, Wn)

(2) 采用汉明窗设计高通FIR滤波器

在b=fir1(n, Wn, 'ftype')中，当ftype=high时，可设计高通滤波器。其语法格式为

b=fir1(n, Wn, 'high')

(3) 采用汉明窗设计带通FIR滤波器

在b=fir1(n, Wn)中，当Wn=[W1  W2]时，fir1函数可得到带通滤波器，其通带为W1 < W < W2

W1 和 W2分别为通带的下限频率和上限频率。其语法格式为

b=fir1(n, [W1  W2])

(4) 采用汉明窗设计带阻FIR滤波器

在b = fir1(n,Wn,'ftype')中，当ftype=stop，Wn=[W1  W2]时，fir1函数可得到带阻滤波器，其语法格式为

b=fir1(n,  [W1  W2], 'stop')

(5) 采用其他窗口函数设计FIR滤波器

使用Window参数，可以用其他窗口函数设计出各种加窗滤波器，Window参数可采用的窗口函数有Boxcar，Hamming，Bartlett，Blackman，Kasier和Chebwin等。其默认时为Hamming窗。例如，采用Bartlett窗设计带阻滤波器，其语法结构为

b=fir1(n, [W1  W2], 'stop', Bartlett[n+1])

    注意：用fir1函数设计高通和带阻滤波器时，所使用的阶数n应为偶数，当输入的阶数n为奇数时，fir1函数会自动将阶数增加1形成偶数。

36.2.2 fir1设计低通滤波器实例

    下面我们通过一个实例来讲解fir1的用法。原始信号是由50Hz正弦波和200Hz的正弦波组成，将200Hz的正弦波当做噪声滤掉，下面通过函数fir1设计一组低通滤波器系数，其阶数是30，截止频率为0.25（也就是125Hz）。Matlab运行代码如下：

1. %****************************************************************************************
   
2. %                             FIR低通滤波器设计
   , H6 J0 l: e( w; Q
3. %***************************************************************************************
   + y. @' e: g  s/ {9 j, s& ^+ w5 t
4. fs=1000;                %设置采样频率 1k
   9 V+ q# _! x9 J( I3 l5 h, N
5. N=1024;                %采样点数      
   & T, t, Z; Q1 w0 o4 [1 _7 \
6. n=0:N-1;
   ( N2 Z" |- i3 Z9 G0 \
7. t=0:1/fs:1-1/fs;          %时间序列
   
8. f=n*fs/N;               %频率序列
   
9. 
   . x' {5 H+ w: u& ~' w# `& F
10. Signal_Original=sin(2*pi*50*t);      %信号50Hz正弦波
    , s. E1 |) H( i/ d
11. Signal_Noise=sin(2*pi*200*t);      %噪声200Hz正弦波
    
12. 
    . ]! D- |- [+ `% }; S! d) ^
13. Mix_Signal=Signal_Original+Signal_Noise;      %将信号Signal_Original和Signal_Original合成一个信号进行采样               
    
14. subplot(221);
    % Z3 I. R! w( l
15. plot(t, Mix_Signal);   %绘制信号Mix_Signal的波形                                                
    $ ~5 m% R% F" i
16. xlabel('时间');
    ) ]4 s2 D; o" ]; N
17. ylabel('幅值');
    . O0 ]/ {1 j/ l4 Y1 o) D+ J
18. title('原始信号');
    2 D8 f; `4 f  @5 u
19. grid on;
    $ g; L4 F' q$ w  F+ w
20.   
    ' Q& Y5 B$ @7 [! G) G4 V
21. subplot(222);
    & Z( t4 e% B4 ]& ~% |) e
22. y=fft(Mix_Signal, N);     %对信号 Mix_Signal做FFT   
    
23. plot(f,abs(y));
      {6 y2 L$ B! V( F
24. xlabel('频率/Hz');
    1 F' n. x8 ^3 |, a; w; i+ y
25. ylabel('振幅');
    4 e/ V; v; T: D
26. title('原始信号FFT');
    & f5 B0 Q! l5 Z2 X+ A9 W
27. grid on;
    
28. 
    1 }& L$ `9 v1 D$ Q3 |# \, U
29. b = fir1(30, 0.25);       %30阶FIR低通滤波器，截止频率125Hz
    
30. %y2= filter(b, 1, x);
    
31. y2=filtfilt(b,1,x);           %经过FIR滤波器后得到的信号
    
32. Ps=sum(Signal_Original.^2);          %信号的总功率
    % z/ ?& [$ E1 H: p( l9 s
33. Pu=sum((y2-Signal_Original).^2);     %剩余噪声的功率
    ; A6 I* p0 ]" _8 `
34. SNR=10*log10(Ps/Pu);               %信噪比
    : B% w" o5 n& X) H( m+ Q* r
35. 
    * R" {- J2 f9 q& @2 e3 b8 P
36. y3=fft(y2, N);            %经过FIR滤波器后得到的信号做FFT
    / J+ r1 G& Z+ n* F- j8 o+ K5 f/ D" V
37. subplot(223);                              
    & X2 C  _8 Z8 s9 u. {+ k( a4 V
38. plot(f,abs(y3));
    
39. xlabel('频率/Hz');
    
40. ylabel('振幅');
    
41. title('滤波后信号FFT');
    
42. grid on;
    
43. 
    
44. [H,F]=freqz(b,1,512);        %通过fir1设计的FIR系统的频率响应
    
45. subplot(224);
    * M8 U* R& x( P) y
46. plot(F/pi,abs(H));           %绘制幅频响应
    
47. xlabel('归一化频率');        
    
48. title(['Order=',int2str(30),'    SNR=',num2str(SNR)]);
    
49. grid on;

复制代码

Matlab的运行结果如下：

从运行结果的FFT和信噪比来看，滤波效果比较明显。

/ B6 i! \6 {' ]0 i$ ^36.2.3 fir1设计高通滤波器实例

    下面我们通过一个实例来讲解fir1的高通滤波器的用法。原始信号是由50Hz正弦波和200Hz的正弦波组成，将50Hz的正弦波当做噪声滤掉，下面通过函数fir1设计一组高通滤波器系数，其阶数是30，截止频率为0.25（也就是125Hz）。Matlab运行代码如下：

1. %****************************************************************************************
   
2. %                             FIR高通滤波器设计
   - `: b0 n2 `$ S+ a/ E8 A
3. %***************************************************************************************
   
4. fs=1000;                 %设置采样频率 1k
   
5. N=1024;                 %采样点数      
   
6. n=0:N-1;
   ( C: g2 z! w1 {
7. t=0:1/fs:1-1/fs;           %时间序列
   
8. f=n*fs/N;                %频率序列
   
9. 
   * ]9 b6 p/ n6 g( r$ }- b
10. Signal_Original=sin(2*pi*200*t);   %信号200Hz正弦波
    , ?* \  @( w& V5 F( m/ C
11. Signal_Noise=sin(2*pi*50*t);      %噪声50Hz正弦波
    
12. 
    & |% ~( U0 I( ]; U7 A0 v
13. Mix_Signal=Signal_Original+Signal_Noise;      %将信号Signal_Original和Signal_Original合成一个信号进行采样               
    
14. subplot(221);
    
15. plot(t, Mix_Signal);   %绘制信号Mix_Signal的波形                                                
    
16. xlabel('时间');
    ; |! g; ~! w4 v/ C- T8 v. W. E
17. ylabel('幅值');
    & n- N/ |. u4 ?) Z" P8 ?; c/ \
18. title('原始信号');
    - t6 ~  w; L' c
19. grid on;
    
20.   
    
21. subplot(222);
    
22. y=fft(Mix_Signal, N);     %对信号 Mix_Signal做FFT   
    
23. plot(f,abs(y));
    
24. xlabel('频率/Hz');
    1 W8 n* @2 q# @% V0 C- U
25. ylabel('振幅');
    
26. title('原始信号FFT');
    ( n) ]6 Y& R) `. S7 ^6 Z
27. grid on;
    
28. 
    
29. b = fir1(30, 0.25, 'high');   %30阶FIR低通滤波器，截止频率125Hz
    2 ^- j/ R* \7 c& w
30. %y2= filter(b, 1, x);
    2 [( _" w8 }( n5 ]# m( y2 K
31. y2=filtfilt(b,1,x);           %经过FIR滤波器后得到的信号
    
32. Ps=sum(Signal_Original.^2);           %信号的总功率
    * {: t, R* o5 I! u7 B5 U
33. Pu=sum((y2-Signal_Original).^2);       %剩余噪声的功率
    ) s; p6 L  s/ P8 L6 h6 A( f! O
34. SNR=10*log10(Ps/Pu);                 %信噪比
    
35. 
    
36. y3=fft(y2, N);            %经过FIR滤波器后得到的信号做FFT
    
37. subplot(223);                              
    
38. plot(f,abs(y3));
    
39. xlabel('频率/Hz');
    
40. ylabel('振幅');
    
41. title('滤波后信号FFT');
    
42. grid on;
    ( `- k& C; L7 u' g  F1 R3 b) e6 U
43. 
    4 K. y8 m9 G5 ^, S% O! I+ T
44. [H,F]=freqz(b,1,512);        %通过fir1设计的FIR系统的频率响应
    6 q4 d1 f% N) W, J- o# y6 ~7 `4 h
45. subplot(224);
    " i) S# D0 l$ p& c
46. plot(F/pi,abs(H));            %绘制幅频响应
    : Q- s8 {0 x  V1 v) m# S' N
47. xlabel('归一化频率');        
    
48. title(['Order=',int2str(30),' SNR=',num2str(SNR)]);
    
49. grid on;

复制代码

Matlab的运行结果如下：

从运行结果的FFT和信噪比来看，滤波效果比较明显。

) X2 u8 L/ d4 u+ y% j% c0 r+ X! O/ `

36.2.4 fir1设计带通滤波器实例

    下面我们通过一个实例来讲解fir1的带通滤波器的用法。原始信号是由50Hz正弦波和200Hz的正弦波组成，设计通带为125Hz到300Hz，下面通过函数fir1设计一组带通滤波器系数，其阶数是30，通带为0.25 < W <0.6。Matlab运行代码如下：

1. %****************************************************************************************
   + J" f8 j' v* T' M( [4 ?' d
2. %                             FIR带通滤波器设计
   
3. %***************************************************************************************
   
4. fs=1000;                 %设置采样频率 1k
   
5. N=1024;                 %采样点数      
   * T/ x6 b* O6 M( Q/ _
6. n=0:N-1;
   , }8 l: i: U" F+ I
7. t=0:1/fs:1-1/fs;           %时间序列
   5 x0 p) [# I, L, ?0 h! [
8. f=n*fs/N;                %频率序列
   6 T( E* ?1 ]3 u: e7 ^4 t2 X, Y5 b0 U8 K) h
9. 
   
10. Signal_Original=sin(2*pi*200*t);     %信号200Hz正弦波
    ) M9 O: @, D6 O2 Y+ |
11. Signal_Noise=sin(2*pi*50*t);        %噪声50Hz正弦波
    
12. 
    
13. Mix_Signal=Signal_Original+Signal_Noise;      %将信号Signal_Original和Signal_Original合成一个信号进行采样               
    4 Y  ^% {3 b( K/ f0 K6 c- |
14. subplot(221);
    / X/ o7 W5 G1 _
15. plot(t, Mix_Signal);   %绘制信号Mix_Signal的波形                                                
    : c5 F( m+ J* P: l- S& m7 r! V
16. xlabel('时间');
    1 C: o! B" u4 ^3 y
17. ylabel('幅值');
    
18. title('原始信号');
    
19. grid on;
    $ a4 L* h: z6 I& b
20.   
    ' l1 h5 W8 h) G' ?6 t* x* m
21. subplot(222);
    
22. y=fft(Mix_Signal, N);     %对信号 Mix_Signal做FFT   
    
23. plot(f,abs(y));
    5 K! P: h  B5 r5 Q% J5 t4 [
24. xlabel('频率/Hz');
    
25. ylabel('振幅');
    ; a8 _6 o5 C' D* \# x) M
26. title('原始信号FFT');
    
27. grid on;
    2 |& e1 U; m: U3 V( S4 w
28. 
    ! o% R/ G3 d! T) M3 O2 Q4 ^
29. b = fir1(30, [0.25 0.6]);   %30阶FIR低通滤波器，截止频率125Hz
    , s" y2 c% p& _
30. %y2= filter(b, 1, x);
    4 z& C0 `7 s6 V- s! |' I4 Z+ v
31. y2=filtfilt(b,1,x);           %经过FIR滤波器后得到的信号
    $ n' [& X- }: M% _' m) _
32. Ps=sum(Signal_Original.^2);           %信号的总功率
    
33. Pu=sum((y2-Signal_Original).^2);       %剩余噪声的功率
    + l. W/ N& Q0 v8 |0 m6 R$ |
34. SNR=10*log10(Ps/Pu);                 %信噪比
    
35. 
    8 C5 b. l% Y( o: h. V/ p. ^
36. y3=fft(y2, N);            %经过FIR滤波器后得到的信号做FFT
    ! l5 P+ q9 l" w1 X2 y+ _9 N3 L1 s
37. subplot(223);                              
    
38. plot(f,abs(y3));
    # j' t$ W* e9 a, ]; i
39. xlabel('频率/Hz');
    9 N/ m' |5 K7 p9 H" D6 b  q" Q
40. ylabel('振幅');
    
41. title('滤波后信号FFT');
    
42. grid on;
    + n/ k. o3 o" X' R4 R
43. 
    3 ]/ O$ C! B7 P4 ^3 \
44. [H,F]=freqz(b,1,512);        %通过fir1设计的FIR系统的频率响应
    % I! v$ v2 X7 |7 F8 X! w
45. subplot(224);
    
46. plot(F/pi,abs(H));            %绘制幅频响应
    
47. xlabel('归一化频率');        
    
48. title(['Order=',int2str(30),' SNR=',num2str(SNR)]);
    
49. grid on;

复制代码

Matlab运行结果如下：

从运行结果的FFT和信噪比来看，滤波效果比较明显。

36.2.5 fir1设计带阻滤波器实例

    下面我们通过一个实例来讲解fir1的带阻滤波器的用法。原始信号是由50Hz正弦波和200Hz的正弦波组成，设计阻带为125Hz到300Hz，下面通过函数fir1设计一组带阻滤波器系数，其阶数是30，阻带为0.25 < W <0.6。Matlab运行代码如下：

1. %****************************************************************************************
   8 t+ T3 a+ T, r) R) t" h. h
2. %                             FIR带阻滤波器设计
   
3. %***************************************************************************************
   , j( T. O7 k' }2 b
4. fs=1000;                 %设置采样频率 1k
   
5. N=1024;                 %采样点数      
   
6. n=0:N-1;
   
7. t=0:1/fs:1-1/fs;           %时间序列
   
8. f=n*fs/N;                %频率序列
   2 g# S# q! o6 j9 m: \
9. 
   / }. a6 D* z3 i$ Z3 s
10. Signal_Original=sin(2*pi*50*t);     %信号50Hz正弦波
    
11. Signal_Noise=sin(2*pi*200*t);      %噪声200Hz正弦波
    7 D" ~  j' s  F4 X! `8 |/ T
12. 
    9 O' ?/ V  c# Z4 h+ h# G
13. Mix_Signal=Signal_Original+Signal_Noise;      %将信号Signal_Original和Signal_Original合成一个信号进行采样               
    / l* V5 W% l7 Z- r2 [
14. subplot(221);
    : o* b% Z" g  Q
15. plot(t, Mix_Signal);   %绘制信号Mix_Signal的波形                                                
    
16. xlabel('时间');
    . @, p7 E$ L1 g! {" P
17. ylabel('幅值');
    ) L0 F# \' `  r0 K, p6 E; F
18. title('原始信号');
    
19. grid on;
    
20.   
    
21. subplot(222);
      p3 e: @, Y& J9 S& H# G7 M0 {8 ]+ F' r
22. y=fft(Mix_Signal, N);     %对信号 Mix_Signal做FFT   
    
23. plot(f,abs(y));
    8 f7 e9 f' n/ \1 V9 s9 X4 a: T
24. xlabel('频率/Hz');
    - o* [" T! ~& U3 t
25. ylabel('振幅');
    
26. title('原始信号FFT');
    , a% c( Q& ]7 [2 b9 {: T
27. grid on;
    
28. 
    
29. b = fir1(30, [0.25 0.6], 'stop');   %30阶FIR低通滤波器，截止频率125Hz
    ' z8 Y0 o/ _. F9 V( D1 k
30. %y2= filter(b, 1, x);
    
31. y2=filtfilt(b,1,x);               %经过FIR滤波器后得到的信号
    9 w5 h6 q) |/ @; v
32. Ps=sum(Signal_Original.^2);   %信号的总功率
    2 ^; r* D# L- @9 \
33. Pu=sum((y2-Signal_Original).^2);     %剩余噪声的功率
    5 S6 v) }! N# a' ]3 k
34. SNR=10*log10(Ps/Pu);               %信噪比
    ; N6 n; B* @" D- T) M1 ^
35. 
    ) J+ N0 Z4 a- S
36. y3=fft(y2, N);            %经过FIR滤波器后得到的信号做FFT
    
37. subplot(223);                              
    / {; m+ z  @+ Z, G7 \4 e
38. plot(f,abs(y3));
    / x7 K2 y" z0 O/ {2 |% u
39. xlabel('频率/Hz');
    
40. ylabel('振幅');
    
41. title('滤波后信号FFT');
    
42. grid on;
    $ O4 ]/ ?2 D+ p1 q
43. 
    
44. [H,F]=freqz(b,1,512);        %通过fir1设计的FIR系统的频率响应
    
45. subplot(224);
      m% t4 @* c: L) m& E2 l; D
46. plot(F/pi,abs(H));            %绘制幅频响应
    9 x: Y& S7 T  L' S6 E- P* @
47. xlabel('归一化频率');        
    $ h% `1 I& V9 H) u4 T$ ]
48. title(['Order=',int2str(30),' SNR=',num2str(SNR)]);
    ' f+ A" }% W9 ]+ c
49. grid on;

复制代码

Matlab运行结果如下：

从运行结果的FFT和信噪比来看，滤波效果比较明显。

. _, r0 E) ~5 K' u$ i) v( j36.2.6 切比雪夫窗口函数设计带通滤波器实例

    下面我们通过一个实例来讲解fir1设计切比雪夫窗口的的带通滤波器。原始信号是由50Hz正弦波和200Hz的正弦波组成，设计通带为125Hz到300Hz，下面通过函数fir1设计一组带通滤波器系数，其阶数是30，通带为0.25 < W <0.6，并且具有25db波纹的切比雪夫窗。Matlab运行代码如下：

1. %****************************************************************************************
   
2. %                             切比雪夫窗口函数设计带通滤波器
   $ P- x5 J" Q$ I1 H8 Y" K7 V
3. %***************************************************************************************
   
4. fs=1000;                 %设置采样频率 1k
   
5. N=1024;                 %采样点数      
   
6. n=0:N-1;
   
7. t=0:1/fs:1-1/fs;           %时间序列
   9 u* {' c( _/ l6 B
8. f=n*fs/N;                %频率序列
   
9. 
   3 v9 ^/ Y+ W+ b1 j! Y6 R5 Y
10. Signal_Original=sin(2*pi*200*t);     %信号200Hz正弦波
    2 h8 a5 i1 h: }3 ~
11. Signal_Noise=sin(2*pi*50*t);      %噪声50Hz正弦波
    ) I$ H, M0 o: ~7 A$ s
12. 
    7 p7 N) t8 h0 c4 _% W. F8 s
13. Mix_Signal=Signal_Original+Signal_Noise;      %将信号Signal_Original和Signal_Original合成一个信号进行采样               
    ' R3 B1 v. a" g, F
14. subplot(221);
    . X  A+ X" R4 U6 t1 ~0 f& b$ y
15. plot(t, Mix_Signal);   %绘制信号Mix_Signal的波形                                                
    
16. xlabel('时间');
    8 c* k5 f! g3 c" `
17. ylabel('幅值');
    
18. title('原始信号');
    
19. grid on;
    
20.   
    & A& k* z0 f) ~: ~4 d
21. subplot(222);
    
22. y=fft(Mix_Signal, N);     %对信号 Mix_Signal做FFT   
    3 S: p9 A1 B) g, x2 K* p* G
23. plot(f,abs(y));
    : I7 d" ~  ^1 s' @3 }7 p/ Y$ A
24. xlabel('频率/Hz');
      w) j1 V; {6 l2 d* S1 X
25. ylabel('振幅');
    
26. title('原始信号FFT');
    8 G0 Y0 m3 |7 E' V$ {) `- Z
27. grid on;
    " W$ t% e1 n8 J8 ?" v. |- p
28. 
    % t" W2 t0 z4 Q# O  ?/ R8 `; A
29. Window = chebwin(31, 25);           %25db的切比雪夫窗
    * A% L7 T% z! H
30. b = fir1(30, [0.25 0.6], Window);   %30阶FIR低通滤波器，截止频率125Hz
    + V( R; w2 }, s& C3 U( W* a
31. %y2= filter(b, 1, x);
    
32. y2=filtfilt(b,1,x);           %经过FIR滤波器后得到的信号
    - s; y  R* ?) a* G! M
33. Ps=sum(Signal_Original.^2);           %信号的总功率
    9 x( w+ R! V( }1 |! y8 N
34. Pu=sum((y2-Signal_Original).^2);       %剩余噪声的功率
    
35. SNR=10*log10(Ps/Pu);                 %信噪比
    1 h& K) \# u$ ]8 }( B, z0 Z
36. 
    
37. y3=fft(y2, N);            %经过FIR滤波器后得到的信号做FFT
    
38. subplot(223);                              
    
39. plot(f,abs(y3));
    
40. xlabel('频率/Hz');
    
41. ylabel('振幅');
    
42. title('滤波后信号FFT');
    9 @: \5 H4 Z  y9 c. m+ C% O! Q
43. grid on;
    % ~- t, J5 \' P: u" O5 A' G  G
44. 
    $ f& u/ ^7 k' o3 P/ N
45. [H,F]=freqz(b,1,512);        %通过fir1设计的FIR系统的频率响应
    * w( z' x( B2 N: T5 ?
46. subplot(224);
    
47. plot(F/pi,abs(H));            %绘制幅频响应
    
48. xlabel('归一化频率');        
    + ?% W9 h! F9 n5 |
49. title(['Order=',int2str(30),' SNR=',num2str(SNR)]);
    - R) G: ]5 I: B$ u: H' X
50. grid on;

复制代码

Matlab运行结果如下：

通过归一化频率可以看出切比雪夫窗口是有一定纹波的。不过从FFT结果和信噪比来看，通过切比雪夫窗口做的滤波效果也是比较明显的。

36.3 fir2函数
36.3.1 fir2函数介绍

    函数fir2用来设计有任意频率响应的各种加窗FIR滤波器。

语法：

    b = fir2(n,f,m)

    b = fir2(n,f,m,window)

    b = fir2(n,f,m,npt)

    b = fir2(n,f,m,npt,window)

    b = fir2(n,f,m,npt,lap)

    b = fir2(n,f,m,npt,lap,window)

参数n为滤波器的阶数。

参数f为频率点矢量，且f[0, 1], f=1对应于0.5fs。矢量f按升序排列，且第一个元素必须为0，最后一个必须为1，并可以包含重复的频率点。

参数m为幅度点矢量，在矢量m中包含了与f相对应的期望得到的滤波器幅度。

参数Window用来指导所使用的窗函数类型，其默认值为汉明窗。

参数npt用来指定fir2函数对频率响应进行内插的点数。

参数lap用来指定fir2函数在重复频率点附近插入的区域大小。

2 e/ u/ y% \8 Q/ j8 M36.3.2 fir2设计低通滤波器

    fir2函数是用来设计任意频率响应的各种加窗FIR滤波器，此函数使用也比较简单，但是要采样的频率点和幅值不好把握，关于这个函数我们仅提供一个低通滤波器的设计。

    原始信号是由50Hz正弦波和200Hz的正弦波组成，将200Hz的正弦波当做噪声滤掉，下面通过函数fir2进行设计。其中频率点矢量和幅度点矢量配置如下：

F = [0 0.1  0.2  0.3  0.4 0.5  0.6  0.7  0.8  1];

A = [1  1  1    1     0 0    0   0    0    0];   

Matlab运行的代码如下：

1. %****************************************************************************************
   
2. %                            fir2设计低通滤波器
   
3. %***************************************************************************************
   
4. fs=1000;                 %设置采样频率 1k
   
5. N=1024;                 %采样点数      
   
6. n=0:N-1;
   
7. t=0:1/fs:1-1/fs;           %时间序列
   
8. f=n*fs/N;                %频率序列
   0 C% `: S6 ?2 ~2 c! a9 B/ {- r/ `
9. 
   : B0 a/ F& P. F! P& c% n6 r* t
10. Signal_Original=sin(2*pi*50*t);     %信号50Hz正弦波
    
11. Signal_Noise=sin(2*pi*200*t);      %噪声200Hz正弦波
    
12. 
    
13. Mix_Signal=Signal_Original+Signal_Noise;      %将信号Signal_Original和Signal_Original合成一个信号进行采样               
    " z/ ~+ U; `/ L5 W8 K
14. subplot(221);
    * x# H$ G) R- X7 t
15. plot(t, Mix_Signal);   %绘制信号Mix_Signal的波形                                                
    
16. xlabel('时间');
    
17. ylabel('幅值');
    
18. title('原始信号');
    ; R0 o# S' N, K# q0 Z
19. grid on;
    
20.   
    
21. subplot(222);
      O9 B5 N: d' W
22. y=fft(Mix_Signal, N);     %对信号 Mix_Signal做FFT   
    
23. plot(f,abs(y));
    4 g$ ]" [5 U0 Q  R
24. xlabel('频率/Hz');
      x" b! T: N+ B
25. ylabel('振幅');
    " A4 ]$ q# \$ U) r, m! x1 r
26. title('原始信号FFT');
    0 H8 e! T3 ?: K
27. grid on;
    
28. 
    
29. F = [0 0.1  0.2  0.3  0.4 0.5  0.6  0.7  0.8  1];  %表示要采样的点
    
30. A = [1  1     1     1     0    0     0     0     0     0];   %表示采样点的幅值
      g1 O. I0 o- F' [
31. b = fir2(30, F, A);   %30阶FIR低通滤波器
    4 ]: k% q- t  R. X& Y4 l. M
32. %y2= filter(b, 1, x);
    ! |4 z& P2 V; r  O5 U0 @9 `
33. y2=filtfilt(b,1,x);           %经过FIR滤波器后得到的信号
    $ R4 \- n$ H9 C) n
34. Ps=sum(Signal_Original.^2);           %信号的总功率
    & l5 @2 D* u, p) `5 E' y) B
35. Pu=sum((y2-Signal_Original).^2);       %剩余噪声的功率
    6 W0 W5 o4 Q* s: z
36. SNR=10*log10(Ps/Pu);                 %信噪比
    
37. 
    2 \5 G# V1 B! H' w$ u+ b
38. y3=fft(y2, N);            %经过FIR滤波器后得到的信号做FFT
    # N, F: c# i7 u! E* A/ g
39. subplot(223);                              
    4 x9 ?1 B# v( F0 L1 s
40. plot(f,abs(y3));
    . T/ i) @5 L/ i2 q& Q+ q6 u2 x
41. xlabel('频率/Hz');
    7 x1 Z$ b/ C* m* N1 E; j
42. ylabel('振幅');
    
43. title('滤波后信号FFT');
    
44. grid on;
    6 }( o! T& s+ `0 |$ I
45. 
    4 \2 S, O( X& ~. n1 j
46. [H,F]=freqz(b,1,512);        %通过fir1设计的FIR系统的频率响应
    " y% s9 v( d3 U( c5 F7 A) z
47. subplot(224);
    3 e% k3 I0 K) O% L+ j" t( `
48. plot(F/pi,abs(H));            %绘制幅频响应
    2 S6 w4 C2 ~6 h0 A, {& k
49. xlabel('归一化频率');        
    
50. title(['Order=',int2str(30),' SNR=',num2str(SNR)]);
    5 }! g% H$ p! D; o8 S
51. grid on;

复制代码

Matlab运行结果如下：

从FFT结果和信噪比来看，fir2任意滤波器设计效果也是比较明显的。

36.4 总结

    本章节主要讲解了函数fir1和函数fir2的使用，想深入的掌握这两个函数，还需要大家多多练习。

4 X2 B( z. E# t" X

支持原创，沙发

感觉要复习一下了

学习学习下

学习。。。。。。

学习一下，，

怎样解决群延时带来的数据丢失呢？