發(fā)布時間：2016-10-10 17:31

  本文關(guān)鍵詞：濾波器，由筆耕文化傳播整理發(fā)布。

基于Matlab的FIR濾波器設(shè)計與實現(xiàn)

一、摘要

　　前面一篇文章介紹了通過FDATool工具箱實現(xiàn)濾波器的設(shè)計，見“基于Matlab中FDATool工具箱的濾波器設(shè)計及相關(guān)文件的生成”，這里通過幾個例子說明采用Matlab語言設(shè)計FIR濾波器的過程。

二、實驗平臺

　　Matlab7.1

三、實驗原理

　　以低通濾波器為例，其常用的設(shè)計指標有：

　　其中，以1、2、3、4條最為常用。5、6條在程序中估算濾波器階數(shù)等參數(shù)時會用到。

　　數(shù)字頻率 = 模擬頻率/采樣頻率

四、實例分析 例1  用凱塞窗設(shè)計一FIR低通濾波器，通帶邊界頻率Ωp=0.3pi，阻帶邊界頻率 Ωs=0.5pi，阻帶衰減δs不小于50dB。 方法一：手動計算濾波器階數(shù)N和β值，之后在通過程序設(shè)計出濾波器。

第一步：通過過渡帶寬度和阻帶衰減，計算濾波器的階數(shù)B和β值。

第二步：通過程序設(shè)計濾波器。

程序如下：

b = fir1(29,0.4,kaiser(30,4.55));

[h1,w1]=freqz(b,1);  

plot(w1/pi,20*log10(abs(h1))); 

axis([0,1,-80,10]); 

grid;

xlabel('歸一化頻率/p') ;

ylabel('幅度/dB') ;

 

波形如下：

 

方法二：

采用[n,Wn,beta,ftype] = kaiserord(f,a,dev)函數(shù)來估計濾波器階數(shù)等，得到凱塞窗濾波器。

　　這里的函數(shù)kaiserord(f,a,dev)或者kaiserord(f,a,dev,fs)：

　　f為對應(yīng)的頻率，fs為采樣頻率；當f用數(shù)字頻率表示時，fs則不需要寫。

　　a=[1 0]為由f指定的各個頻帶上的幅值向量，一般只有0和1表示；a和f長度關(guān)系為（2*a的長度）- 2=（f的長度）

　　devs=[0.05 10^(-2.5)]用于指定各個頻帶輸出濾波器的頻率響應(yīng)與其期望幅值之間的最大輸出誤差或偏差，長度與a相等，計算公式：

阻帶衰減誤差=αs，通帶衰減誤差=αp，可有濾波器指標中的3、4條得到。

　　fs缺省為2Hz。

 

程序如下：

fcuts = [0.3  0.5]; %歸一化頻率omega/pi，這里指通帶截止頻率、阻帶起始頻率

mags = [1 0];

devs = [0.05 10^(-2.5)];

[n,Wn,beta,ftype] = kaiserord(fcuts,mags,devs);  %計算出凱塞窗N，，beta的值

hh = fir1(n,Wn,ftype,kaiser(n+1,beta),'noscale'); 

freqz(hh);

 

波形如下：

 

　　實際中，一般調(diào)用MATLAB信號處理工具箱函數(shù)remezord來計算等波紋濾波器階數(shù)N和加權(quán)函數(shù)W(ω)，調(diào)用函數(shù)remez可進行等波紋濾波器的設(shè)計，直接求出濾波器系數(shù)。函數(shù)remezord中的數(shù)組fedge為通帶和阻帶邊界頻率，數(shù)組mval是兩個邊界處的幅值，而數(shù)組dev是通帶和阻帶的波動，fs是采樣頻率單位為Hz。

 

例2  利用雷米茲交替算法設(shè)計等波紋濾波器，設(shè)計一個線性相位低通FIR數(shù)字濾波器，其指標為：通帶邊界頻率fc=800Hz，阻帶邊界fr=1000Hz，通帶波動 阻帶最小衰減At=40dB，采樣頻率fs=4000Hz。 

 

解：在MATLAB中可以用remezord 和remez兩個函數(shù)設(shè)計

 

程序如下： 

 

fedge=[800 1000]; 

 

mval=[1 0]; 

 

dev=[0.0559 0.01];

 

fs=4000; 

 

[N,fpts,mag,wt]=remezord(fedge,mval,dev,fs); 

 

b=remez(N,fpts,mag,wt); 

 

[h,w]=freqz(b,1,256); 

 

plot(w*2000/pi,20*log10(abs(h))); 

 

grid; 

 

xlabel('頻率/Hz') ;

 

ylabel('幅度/dB');

 

波形如下：

例3 利用MATLAB編程設(shè)計一個數(shù)字帶通濾波器，指標要求如下：通帶邊緣頻率：Ωp1=0.45pi，Ωp2=0.65pi，通帶峰值起伏：δ1<=1[dB]。阻帶邊緣頻率：Ωs1=0.3pi，Ωs2=0.8pi，最小阻帶衰減：δ2>=40[dB] 。 方法一：窗函數(shù)法

 

程序如下：

[n,wn,bta,ftype]=kaiserord([0.3 0.45 0.65 0.8],[0 1 0],[0.01 0.1087 0.01]);%用kaiserord函數(shù)估計出濾波器階數(shù)n和beta參數(shù)

h1=fir1(n,wn,ftype,kaiser(n+1,bta),'noscale');

[hh1,w1]=freqz(h1,1,256);

figure(1)

subplot(2,1,1)

plot(w1/pi,20*log10(abs(hh1)))

grid

xlabel('歸一化頻率w');ylabel('幅度/db');

subplot(2,1,2)

plot(w1/pi,angle(hh1))

grid

xlabel('歸一化頻率w');ylabel('相位/rad');

 

波形如下：

 

濾波器系數(shù)為：

h1 =

 

  Columns 1 through 8 

 

    0.0041    0.0055   -0.0091   -0.0018   -0.0056   -0.0000    0.0391   -0.0152

 

  Columns 9 through 16 

 

   -0.0381    0.0077   -0.0293    0.0940    0.0907   -0.2630   -0.0517    0.3500

 

  Columns 17 through 24 

 

   -0.0517   -0.2630    0.0907    0.0940   -0.0293    0.0077   -0.0381   -0.0152

 

  Columns 25 through 31 

 

    0.0391   -0.0000   -0.0056   -0.0018   -0.0091    0.0055    0.0041

 

如果直接用freqz(h1,1,256)，得幅頻特性和相頻特性曲線：

  方法二：等波紋法設(shè)計

 

程序如下：

[n,fpts,mag,wt]=remezord([0.3 0.45 0.65 0.8],[0 1 0],[0.01 0.1087 0.01]);%用remezord函數(shù)估算出remez函數(shù)要用到的階n、歸一化頻帶邊緣矢量fpts、頻帶內(nèi)幅值響應(yīng)矢量mag及加權(quán)矢量w，使remez函數(shù)設(shè)計出的濾波器滿足f、a及dev指定的性能要求。
h2=remez(n,fpts,mag,wt);%設(shè)計出等波紋濾波器
[hh2,w2]=freqz(h2,1,256);
figure(2)
subplot(2,1,1)
plot(w2/pi,20*log10(abs(hh2)))
grid
xlabel('歸一化頻率w');ylabel('幅度/db');
subplot(2,1,2)
plot(w2/pi,angle(hh2))
grid
xlabel('歸一化頻率w');ylabel('相位/rad');
h2

 

波形如下：

濾波器系數(shù)如下：

h2 =

 

  Columns 1 through 9 

 

   -0.0013    0.0092   -0.0255   -0.0642    0.1177    0.0922   -0.2466   -0.0466    0.3116

 

  Columns 10 through 17 

 

   -0.0466   -0.2466    0.0922    0.1177   -0.0642   -0.0255    0.0092   -0.0013

 

如果直接用freqz(h2,1,256);得幅頻特性和相頻特性曲線：

 

方法三：采用FDATool工具

 這種方法需要事先計算出濾波器的階數(shù)，bate值，然后設(shè)置相應(yīng)參數(shù)，最后生成濾波器。

 設(shè)置界面如下圖所示：

 

 　　將上述圈圈的區(qū)域設(shè)置好之后，生成濾波器，最后通過analysis菜單可以觀察生成的濾波器的各種特性曲線和濾波器系數(shù)。這里的濾波器系數(shù)跟方法一的一樣。

 

波形如下：

 

五、結(jié)果分析 5.1  濾波器設(shè)計總結(jié)

　　FIR濾波器實現(xiàn)一般采用窗函數(shù)法和等紋波設(shè)計法。窗函數(shù)法還包含兩個分支，一種是用公式先手動算出N值和其他對應(yīng)得窗函數(shù)參數(shù)值，再代入窗函數(shù)和fir1實現(xiàn)，一種是用函數(shù)*rord估算出N和相應(yīng)參數(shù)再用fir1實現(xiàn)。不過要注意*rord會低估或高估階次n，可能會使濾波器達不到指定的性能，這時應(yīng)稍微增加或降低階次。如果截止頻率在0或Nyquist頻率附近，或者設(shè)定的dev值較大，則得不到正確結(jié)果。

　　濾波器實現(xiàn)形式及特點：由于一般的濾波器在利用窗函數(shù)是其通帶波紋和阻帶波紋不同（一般為第一個阻帶波紋最大）因此，在滿足第一個阻帶衰減旁瓣時，比其頻率高的旁瓣，它們的衰減都大大超出要求。而根據(jù)阻帶衰減與項數(shù)的近似關(guān)系N = P(δ2)*fs/TW，可得當阻帶衰減越大，所需項數(shù)越多。

5.2  窗函數(shù)法和等波紋設(shè)計的不同之處

　　窗函數(shù)設(shè)計是通過最小平方積分辦法來設(shè)計的，即該濾波器的誤差為：

　　

　　即要求最小方法來設(shè)計濾波器，這樣的濾波器更忠實于理想濾波器（即濾波系數(shù)更接近于理想濾波器）。

證明如下：

　　

因此，幅度頻譜差值越小，實際濾波器就越接近理想濾波器。

　　而等波紋濾波器是通過最大加權(quán)誤差最小化來實現(xiàn)，其誤差為：

　　

要求該誤差最小來實現(xiàn)濾波器，得出來的濾波系數(shù)較窗函數(shù)設(shè)計相差較遠。

以下通過對例3中的h1及h2作比較。

%sigsum是用來對數(shù)組各元素進行求和

function y=sigsum(n1,n2,n,x);

y=0; 

for i=n1+1-min(n):n2+1-min(n)

    y=y+x(i);

end

 

n=0.001:30.001;

h=2*cos(0.55*pi*(n-15)).*sin(0.175*pi*(n-15))./(pi*(n-15));

delta1=h-h1;

n=0.001:16.001;

h=2*cos(0.55*pi*(n-15)).*sin(0.175*pi*(n-15))./(pi*(n-15));

delta2=h-h2;

y1=sigsum(0,30,[0:30],(abs(delta1).^2))/31;

y2=sigsum(0,16,[0:16],(abs(delta2).^2))/17;

 

 結(jié)果如下：

y1 =

 

  1.9099e-004

 

y2 =

 

    0.0278

 

　　由此得到用窗函數(shù)實現(xiàn)的濾波系數(shù)比用等波紋濾波器系數(shù)的每一項更接近于理想濾波器（y1為用窗函數(shù)實現(xiàn)的與理想濾波器的差值，y2為用等波紋濾波器實現(xiàn)的與理想濾波器的差值）；

 

 

　　對比二者的幅度頻譜可知，等波紋濾波器阻帶邊緣比用窗函數(shù)實現(xiàn)的更平滑（理想濾波器為垂直下降的）。

　　從設(shè)計的角度考慮，由于窗函數(shù)設(shè)計法都是通過已有的窗函數(shù)對理想濾波器的改造，因此，可以用手算的辦法方便的設(shè)計濾波器。

而等波紋濾波器，其實現(xiàn)是通過大量的迭代運算來實現(xiàn)，這樣的方法一般只能通過軟件來設(shè)計。

　　項數(shù)的問題由于等波紋濾波器能較平均的分布誤差，因此對于相同的阻帶衰減，其所需的濾波系數(shù)比窗函數(shù)的要少。

5.3 幾點說明

1.相頻特性曲線形狀不同說明

 

2.調(diào)用firl或者reme函數(shù)時，用scale（缺省方式）對濾波器進行歸一化，即濾波器通帶中心頻率處的響應(yīng)幅值為0db。用noscale不對濾波器歸一化。

 

 

 

 

 

 

 

posted @

  本文關(guān)鍵詞：濾波器，由筆耕文化傳播整理發(fā)布。

本文編號：136465