數字信號處理實驗（六）——FIR濾波器的設計

一、四種線性相位FIR濾波器的振幅響應

1、自編函數

[Hr,w,a,L]=hr_type1(h)（P256）  % h偶對稱，N為奇數，h(n)=h(N-1-n)
[Hr,w,a,L]=hr_type2(h) （P257） % h偶對稱，N為偶數，h(n)=h(N-1-n)
[Hr,w,a,L]=hr_type3(h) （P257） % h奇對稱，N為奇數,h(n)=-h(N-1-n)
[Hr,w,a,L]=hr_type4(h) （P257） % h奇對稱，N為偶數,h(n)=-h(N-1-n)

2、一個demo

clear all;close all;clc
addpath(genpath(pwd)); % 添加當前文件夾下所有路徑

%% 
% 題1． 線性相位 FIR 濾波器的特性：
%% (2)已知濾波器的系統函數如下所示，用以上已編好的函數，確定濾波器的振
% 幅響應Hr(w)以及零點位置：
h1_n = [-4,1,-1,-2,5,6,5,-2,-1,1,-4];     % 偶對稱,N=11
h2_n = [-4,1,-1,-2,5,6,6,5,-2,-1,1,-4];   % 偶對稱,N=12
h3_n = [-4,1,-1,-2,5,0,-5,2,1,-1,4];      % 奇對稱,N=11
h4_n = [ -4,1,-1,-2,5,6,-6,-5,2,1,-1,4];  % 奇對稱,N=12

new_figure('線性相位FIR濾波器');
subplot(2,2,1);
[Hr,w]=hr_type1(h1_n);
plot(w/pi, abs(Hr));
title('第一類線性相位濾波器');
xlabel('w/pi');

subplot(2,2,2);
[Hr,w]=hr_type2(h2_n);
plot(w/pi, abs(Hr));
title('第二類線性相位濾波器');
xlabel('w/pi');

subplot(2,2,3);
[Hr,w]=hr_type3(h3_n);
plot(w/pi, abs(Hr));
title('第三類線性相位濾波器');
xlabel('w/pi');

subplot(2,2,4);
[Hr,w]=hr_type4(h4_n);
plot(w/pi, abs(Hr));
title('第四類線性相位濾波器');
xlabel('w/pi');

 

二、窗函數法

1、窗函數設計參考指標

2、窗函數設計方法一：

（1）根據實際阻帶衰減指標，來確定所使用的窗函數。Matlab提供了幾個函數來實現這些窗函數。

   W=boxcar(N);            % 矩形窗         －21dB
   W=triang(N);            % 三角窗         －25dB
   W=hanning(N);           % 漢寧窗         －44dB
   W=hamming(N);           % 海明窗         －53dB
   W=blackman(N);          % 布萊克曼窗     －74dB
   W=kaiser(N,beta);       % 凱澤窗         －80dB

（2）根據過渡帶來計算出N值。

例： 通帶截止頻率wp, 阻帶截止頻率ws，已知為漢寧窗：N=3.1*2*pi/(ws-wp);

 

（3）求出理想的hd(n)。可使用ideal_lp（P185）來實現理想低通濾波器的沖激響應。

例：由ideal_lp來實現理想帶阻濾波器的沖激響應。

hd=ideal_lp(wc1,N)+ideal_lp(pi,N)-ideal_lp(wc2,N);

 

（4）求得所設計的FIR濾波器的單位抽樣響應:

h=hd.*w

 

（5）一個demo：

clear all;
wp=0.2*pi;
Rp=0.25;
ws=0.3*pi;
As=50;

wd = ws-wp;

N = ceil(6.6*pi/wd);
wn = (wp+ws)/2;

%
hd = ideal_lp(wn,N+1);
w_ham=(hamming(N+1))';
h = hd.*w_ham;
figure(2)
freqz(h,1,512)

 

 

3、窗函數設計方法二：

（1）根據實際阻帶衰減指標，來確定所使用的窗函數。

（2）根據過渡帶來計算出N值。

（3）利用Matlab所提供的函數fir1，來實現fir濾波器。

h=fir1(N,wn,’ftype’,windows(N+1))

›對於高通濾波器和帶阻濾波器，N必須為偶數，N+1為奇數

›‘ftype’指的是：’low’,’bandpass’,’high’,’stop’

›設計的濾波器為N階

（4）一個demo

%example6.2.3
clear all;
wp=0.2*pi;
Rp=0.25;
ws=0.3*pi;
As=50;

wd = ws-wp;
N = ceil(6.6*pi/wd);
wn = (wp+ws)/2;
b=fir1(N,wn/pi,hamming(N+1));
figure(1)
freqz(b,1,512)

 

三、頻率抽樣法：

1、根據給定的N值和過渡點來求得Hd(k) 。

2、利用DFT反變換求的所需的h(n)。

h=real(ifft(H,N));

 

clear all;
N=40;
T1=0.5925;
T2=0.1099;

alpha=(N-1)/2;
l=0:N-1;
wl=(2*pi/N)*l;

hrs=[zeros(1,5),T2,T1,ones(1,7),T1,T2,zeros(1,9),T2,T1,ones(1,7),T1,T2,zeros(1,4)];
k1 = 0:floor((N/2-1));
k2 = floor(N/2)+1:N-1;
angh = [-alpha*(2*pi)/N*k1,0,alpha*(2*pi)/N*(N-k2)];
H=hrs.*exp(j*angh);
h=real(ifft(H,N));

figure()
freqz(h,1,512)
[db,mag,pha,grd,w]=freqz_m(h,1);
[hr,ww,a,L]=hr_type2(h);

figure()
subplot(2,2,1),stem(wl(1:21)/pi,hrs(1:21)),title('H(k)')
subplot(2,2,2),stem(l,h),title('沖激響應')
subplot(2,2,3),plot(ww/pi,hr);grid;title('幅度響應')
subplot(2,2,4),plot(w/pi,db),axis([0,1,-80,10]),grid,title('幅頻響應(db)')