FIR
濾波器的設(shè)計(jì)方法主要有三種:窗函數(shù)法、頻率取樣法、切比雪夫等波紋逼近法。而在這個(gè)實(shí)驗(yàn)里,主要是采取第一種方法。
我們要設(shè)計(jì)一個(gè)濾波器,跟前面設(shè)計(jì)IIR濾波器一樣,得先知道一些關(guān)于濾波器的指標(biāo)。在用窗函數(shù)設(shè)計(jì)FIR濾波器需要知道的指標(biāo)是:通帶,阻帶的邊界頻率,阻帶衰減和通帶衰減。
因?yàn)槭荈IR,有限沖激響應(yīng)的濾波器,我們通常的理想濾波器的單位脈沖響應(yīng) h 是無限長(zhǎng)的,所以我們需要通過窗來截?cái)嗨?,然后?duì)它進(jìn)行濾波器設(shè)計(jì)。
一、 加窗方法設(shè)計(jì)的步驟大概分以下幾步:
1、根據(jù)阻帶的衰減,選擇合適的窗:
最小阻帶衰減 過渡帶帶寬△w
矩形窗 20.9dB 0.92π/M
漢寧窗 43.9dB 3.11π/M
海明窗 54.5dB 3.32π/M
布萊克曼窗 75.3dB 5.56π/M
不同的窗有不同的性質(zhì):不同的窗函數(shù),產(chǎn)生泄漏的大小不一樣,頻率分辨能力也不一樣。信號(hào)的截?cái)喈a(chǎn)生了能量泄漏,而用FFT算法計(jì)算頻譜又產(chǎn)生了柵欄效應(yīng),從原理上講這兩種誤差都是不能消除的,但是我們可以通過選擇不同的窗函數(shù)對(duì)它們的影響進(jìn)行抑制。(矩形窗主瓣窄,旁瓣大,頻率識(shí)別精度最高,幅值識(shí)別精度最低;布萊克曼窗主瓣寬,旁瓣小,頻率識(shí)別精度最低,但幅值識(shí)別精度最高)。
2、根據(jù)窗函數(shù)得到的序列經(jīng)過firl或fir2得到一個(gè)濾波器傳輸函數(shù)系數(shù)的序列。
1)fir1 : 用來設(shè)計(jì)傳統(tǒng)的低通,高通,帶通,帶阻,多頻帶FIR濾波器;
調(diào)用格式:b = fir1(N,Wn);
b = fir1(N,Wn,‘high’);
b = fir1(N,Wn, ‘stop’);
參數(shù)說明:N:階次,濾波器長(zhǎng)度為N+1;
Wn:通帶截止頻率,其值在0~1之間,1對(duì)應(yīng)Fs/2;
b:濾波器系數(shù)。
在上述所有格式中,若不指定窗函數(shù)的類型,fir1自動(dòng)選擇Hamming窗。
2)fir2 : 用來設(shè)計(jì)具有任意幅度響應(yīng)的FIR濾波器。
調(diào)用格式:b = fir2(N, F, M);
參數(shù)說明:F是頻率向量,其值在0~1之間;
M是和F相對(duì)應(yīng)的所希望的幅頻相應(yīng)。
如同fir1, 缺省時(shí)自動(dòng)選用Hamming窗。
3)為了觀測(cè)到設(shè)計(jì)出來的濾波器的特性,我們可以用impz得到它的脈沖響應(yīng),用freqz得到頻率響應(yīng)。
其中在畫頻率響應(yīng)的時(shí)候我們分為幅度和相位畫出。
又因?yàn)槲覀円^測(cè)的是衰減的大小程度,以dB為單位,所以我們?cè)诋嫹鹊臅r(shí)候縱坐標(biāo)應(yīng)該轉(zhuǎn)換成dB。
而相位,由于計(jì)算機(jī)中反正切函數(shù)規(guī)定,在一、二象限中的角度為0到π,三、四象限的角度為0到-π。我們一般要反應(yīng)的角度變化是0到2π,但如果不經(jīng)過處理,實(shí)際結(jié)果會(huì)在π處發(fā)生跳變,跳變的幅度為2π,這就是相位的卷繞。所以我們用解卷繞函數(shù)unwarp(w),使相位在π處不發(fā)生跳變,從而反應(yīng)出真實(shí)的相位變化。
二、實(shí)例
%%第一題、通帶截止頻率為0.2pi,阻帶截止頻率為0.3pi,阻帶衰減不小于50dB,通帶衰減不大于3dB,設(shè)計(jì)一個(gè)FIR濾波器,并驗(yàn)證。
wp=0.2*pi; ws=0.3*pi; %性能指標(biāo)
wdelta=ws-wp; %過渡帶寬度
M=ceil(3.32*pi/wdelta); %濾波器長(zhǎng)度,朝正無窮方向舍入
N=2*M+1; %窗口長(zhǎng)度
wc=(ws+wp)/2; %截止頻率
win=hamming(N); %因?yàn)樗p不小于50dB,所以選擇海明窗,這里得到海明窗的時(shí)域響應(yīng)
b=fir1(N-1,wc/pi,win);
n=0:1:N;
[hi t]=impz(b,1,n);%得到脈沖響應(yīng)
[hf w]=freqz(b,1,512); %得到頻率響應(yīng)
figure(1);
subplot(3,1,1); stem(n,hi);
xlabel('n'); ylabel('Magnitude'); title('impulse response');
subplot(3,1,2); plot(w/pi,20*log10(abs(hf)));
xlabel('Frequency(Hz)'); ylabel('Magnitude(dB)');
title('Frequency response');
subplot(3,1,3); plot(w/pi,180/pi*unwrap(angle(hf)));
xlabel('Frequency(Hz)'); ylabel('Phase(degrees)');
title('Frequency response');
%驗(yàn)證:
nn=0:50;
x1=sin(nn*0.2*pi); x2=sin(nn*0.8*pi);%假設(shè)兩個(gè)信號(hào),低頻和高頻
in=x1+x2; out=filter(b,1,in); %濾波過程
figure(2);
subplot(2,2,1); stem(x1);
xlabel('n'); ylabel('Magnitude');
title('x1'); axis([0 50 -1 1]);
subplot(2,2,2); stem(x2);
xlabel('n'); ylabel('Magnitude');
title('x2'); axis([0 50 -1 1]);
subplot(2,2,3); stem(in);
xlabel('n'); ylabel('Magnitude');
title('x1+x2 before filter'); axis([0 50 -2 2]);
subplot(2,2,4); stem(out);
xlabel('n'); ylabel('Magnitude');
title('x1+x2 after filter'); axis([0 50 -1 1]);
在這題要主要的是求分貝的公式:20*log10(abs(hf)),還有解卷繞函數(shù)的調(diào)用180/pi*unwrap(angle(hf))。
從上面的程序可以得到圖:
%第二題,取1中相同的濾波器類型和階數(shù),通帶和阻帶的臨界頻率也相同,分別用矩形窗,漢寧窗,海明窗和布萊克曼窗設(shè)計(jì)濾波器,并分析它們的性能。
wp=0.2*pi;
ws=0.3*pi; %性能指標(biāo)
wdelta=ws-wp; %過渡帶寬度
M1=ceil(0.92*pi/wdelta);
M2=ceil(3.11*pi/wdelta);
M3=ceil(3.32*pi/wdelta);
M4=ceil(5.56*pi/wdelta);
N1=2*M1+1;
N2=2*M2+1;
N3=2*M3+1;
N4=2*M4+1;
wc=(ws+wp)/2; %截止頻率
win1=boxcar(N1); %矩形窗的時(shí)域響應(yīng)
win2=hanning(N2); %漢寧窗的時(shí)域響應(yīng)
win3=hamming(N3); %海明窗的時(shí)域響應(yīng)
win4=blackman(N4); %布萊克窗的時(shí)域響應(yīng)
b1=fir1(N1-1,wc/pi,win1);
b2=fir1(N2-1,wc/pi,win2);
b3=fir1(N3-1,wc/pi,win3);
b4=fir1(N4-1,wc/pi,win4);
[h1 w1]=freqz(b1,1,512);
[h2 w2]=freqz(b2,1,512) ;
[h3 w3]=freqz(b3,1,512) ;
[h4 w4]=freqz(b4,1,512) ;
figure(3);
subplot(4,2,1); plot(w1/pi,20*log10(abs(h1)));
ylabel('Magnitude'); title('Rectangle Window');
subplot(4,2,2); plot(w1/pi,180/pi*unwrap(angle(h1)));
ylabel('Phase(degrees)'); title('Rectangle Window');
subplot(4,2,3); plot(w2/pi,20*log10(abs(h2)));
ylabel('Magnitude'); title('Hanning Window');
subplot(4,2,4); plot(w2/pi,180/pi*unwrap(angle(h2)));
ylabel('Phase'); title('Hanning Window');
subplot(4,2,5); plot(w3/pi,20*log10(abs(h3)));
ylabel('Magnitude'); title('Hamming Window');
subplot(4,2,6); plot(w3/pi,180/pi*unwrap(angle(h3)));
ylabel('Phase'); title('Hamming Window');
subplot(4,2,7); plot(w4/pi,20*log10(abs(h4)));
ylabel('Magnitude'); title('Blankman Window');
subplot(4,2,8); plot(w4/pi,180/pi*unwrap(angle(h4)));
ylabel('Phase'); title('Blankman Window');
這題只是單純的重復(fù)4次濾波器的設(shè)計(jì),設(shè)計(jì)的方法同第一題一樣。得到的圖如下:
從上可以觀察到,
1、衰減:從上到下,衰減依次增大,以布萊克曼窗的為最優(yōu),矩形窗最差。
2、過渡帶:從上到下,過渡帶逐漸增寬,以矩形窗的為最優(yōu),布萊克曼窗最差。
其實(shí),改善阻帶衰減的一種辦法是加寬過渡帶寬,以犧牲過渡帶換取阻帶衰減的增加。也就是以增加主瓣寬度為代價(jià)來降低旁瓣。
所以,我們?cè)谠O(shè)計(jì)濾波器的時(shí)候,要根據(jù)不同的指標(biāo),性能,需求去選擇合適的窗來進(jìn)行設(shè)計(jì)。