ĐẠI HỌC ĐÀ NẴNG
oOo
TIỂU LUẬN MÔN HỌC
XỬ LÝ TÍN HIỆU SỐ
ĐỀ TÀI:
THIẾT KẾ BỘ LỌC FIR
BẰNG PHƯƠNG PHÁP LẤY MẪU TẦN SỐ
GVHD : TS. NGÔ VĂN SỸ
NHÓM 2 : LÊ PHƯỚC SINH
: PHẠM TUẤN SƠN
: HUỲNH VĂN THỦY
: VŨ BẢO TOÀN
: VÕ ANH TUẤN
LỚP : TĐH K24
Đà Nẵng, tháng 4/2012
DSP – Thiết kế bộ lọc FIR bằng phương pháp lấy mẫu tần số
LỜI MỞ ĐẦU
Xử lý tín hiệu số (Digital Signal Processing – DSP) là một lĩnh vực mới, đã phát
triển mạnh mẽ về mặt lý thuyết và công nghệ. Ngày nay, nó đã trở thành môn học
không thể thiếu cho nhiều ngành học về khoa học và kỹ thuật. Các thuật toán của nó
đã mang lại những thuận tiện cho công nghệ và ứng dụng mới, nhất là các lĩnh vực về
kỹ thuật điện tử và điều khiển tự động.
Phép xử lý cơ bản nhất của DSP là lọc, và các hệ thống được đề cập đến nhiều nhất
trong xử lý tín hiệu số là các bộ lọc số (Digital Filter).
Nếu xét về đáp ứng xung có thể chia các bộ lọc số thành 2 loại chính là bộ lọc có
đáp ứng xung hữu hạn FIR (Finite Impulse Response) còn gọi là lọc không đệ quy, và
bộ lọc có đáp ứng xung vô hạn IIR (Infinte Impulse Response) còn gọi là lọc đệ quy.
Xét về đáp ứng tần số biên độ có thể chia các bộ lọc, FIR hay IIR, thành 4 loại cơ
bản: thông thấp, thông cao, thông dải và chắn dải. Các bộ lọc này có thể được thiết kế
bằng những phương pháp sau đây: Phương pháp cửa sổ (Window Design Techniques),
Phương pháp lấy mẫu tần số (Frequency Sampling Design Techniques) và Phương

Ở trong nhiều ứng dụng như xử lý tiếng nói hoặc xử lý âm thanh, bộ lọc số được
dùng để thực hiện các thao tác chọn tần. Do đó, các chỉ tiêu được đòi hỏi trong vùng
tần số ở mỗi chu kỳ về số lượng và đáp ứng pha của bộ lọc mong muốn. Nói chung
đáp ứng pha cần tuyến tính trong dải thông.
- Trong trường hợp các bộ lọc FIR , Có thể đạt được chính xác yêu cầu về pha tuyến
tính.
- Trong trường hợp các bộ lọc IIR, một dải thông có pha tuyến tính là rất khó đạt.
Do đó, chúng ta chỉ xét các chỉ tiêu về biên độ.
Có 2 nhóm chỉ tiêu:
•
Các chỉ tiêu tuyệt đối (Absolute Specifications) : Đưa ra một tập các yêu cầu
trên hàm đáp ứng biên độ |H(e
jw
)|. Những chỉ tiêu này được sử dụng chung
cho các bộ lọc FIR.
•
Các chỉ tiêu tương đối (Relative Specifications - DB) : Đưa ra các yêu cầu tính
theo decibels (dB), được cho bởi :
Phương pháp này được là một phương pháp phổ biến nhất trong thực tế và được sử
dụng cho cả bộ lọc FIR và IIR.
Trong những
phần tiếp sau đây, để mô tả các chỉ tiêu của bộ lọc chúng ta sẽ xem xét việc thiết kế
một bộ lọc thông thấp như một ví dụ cơ sở để thiết cho bộ lọc FIR.
1.2. Các chỉ tiêu
Trên hình (1.1) là mô tả các chỉ tiêu của bộ lọc FIR thông thấp (Low Pass Filter):
1.2.1.Các chỉ tiêu tuyệt đối
•
Band [0, w
p
] được gọi là dải thông, và δ

] được gọi là dải chuyển tiếp, và không có ràng buộc nào về đáp ứng biên
độ trong dải này.
1.2.2. Các chỉ tiêu tương đối (DB)
R
p
: Độ gợn sóng trong dải thông tính theo dB.
A
s
: Suy hao trong dải chắn tính theo dB.
Quan hệ giữa các chỉ tiêu này như sau:

Tại sao tập trung vào bộ lọc thông thấp?
Các chỉ tiêu trên được đưa ra đối với bộ lọc FIR thông thấp, và tất nhiên đối với các
bộ lọc khác như thông cao HPF (High Pass Filter), thông dải BPF (Band Pass Filter)
đều có thể được định nghĩa tương tự. Tuy nhiên, các tham số thiết kế quan trọng nhất
là các dung sai dải tần và các tần số cạnh-dải. Bởi vậy, trong phần 1 về cơ sở lý thuyết
này chúng ta chỉ tập trung vào bộ lọc FIR thông thấp. Việc thiết kế cụ thể cho bộ lọc
FIR thông dải bằng kỹ thuật cửa sổ sẽ được phát triển trên cơ sở lọc thông thấp và sẽ
được mô tả chi tiết trong phần 2.
Các thuận lợi trong thiết kế và thi hành lọc số FIR
4
(1.2)
(1.3)
0
1
1
2
10
log20 >
+


Độ gợn dải thông
Độ gợn dải chắn
Dải chuyển tiếp
δ
2

R
p

A
s

w
p

w
s

π
0
0
|H(e
jw
)|
Decibels
Hình (1.1) Các chỉ tiêu của bộ lọc FIR
DSP – Thiết kế bộ lọc FIR bằng phương pháp lấy mẫu tần số
Việc thiết kế và thực hiện lọc FIR có những thuận lợi sau đây:
•

−
−
=+++=
1M
0n
n
n
M1
1M
1
10
zbzbzbb)z(H *
(1.4)
Như vậy đáp ứng xung h(n) là:



−≤≤
=
else
Mnb
nh
n
0
10
)(
(1.5)
Và phương trình sai phân là:
)1()1()()(
110






+++=+++=
−
−
−−
−
− M1
0
1M
1
0
1
0
M1
1M
1
10
z
b
b
z
b
b
1bzbzbb)z(H **
(1.8)


t heo tần số, nghĩa là:
αω−β=∠
ω
)e(H
j
π≤ω≤π−
(1.9)
Trong đó :
0=β
hoặc
2
π
±
và
α
là một hằng số.
Đối với bộ lọc FIR nhân quả có đáp ứng xung trong khoảng [0, M-1], thì các điều
kiện tuyến tính là:
1Mn0,0);n1M(h)n(h −≤≤=β−−=
(1.10)
1Mn0,2/);n1M(h)n(h −≤≤π±=β−−−=
(1.11)
Xét phương trình sai phân được cho trong phương trình (1.6) với đáp ứng xung đối
xứng trong phương trình (1.10), ta có:
)1Mn(xb)2Mn(xb)1n(xb)n(xb)n(y
0110
+−++−++−+= *

*++−+−++−+= )]2Mn(x)1n(x[b)]1Mn(x)n(x[b
10

π≤ω<π−=
∑
−
=
ω−ω
,e)n(h)e(H
1M
0n
njj
(1.13)
1.4.1. Đáp ứng xung h(n)
Chúng ta có thể đưa ra ràng buộc pha tuyến tính:
π≤ω<π−αω−=∠
ω
,)e(H
j
(1.14)
trong đó: α là một hằng số trễ pha. Ta đã biết rằng h(n) phải đối xứng, nghĩa là:
2
1M
,1Mn0),n1M(h)n(h
−
=α−≤≤−−=
(1.15)
Do đó h(n) là đối xứng theo α, là chỉ số đối xứng. Có hai kiểu đối xứng:
7

b
0
z


b
1
z
-
1

b
2
b
3
y(n)
x(n)
z
-
1

z
-
1

z
-
1

z
-
1

M=7

π≤ω<π−
(1.16)
Đáp ứng pha là đường thẳng nhưng không đi qua gốc. Trong trường hợp này α
không phải là hằng số trễ pha, nhưng:
α−=
ω
∠
ω
d
)e(Hd
j
(1.17)
là hằng số, chính là trễ nhóm (α là một hằng số trễ nhóm). Trong trường hợp này, các
tần số được làm trễ với một tốc độ không đổi.
Đối với kiểu pha tuyến tính này, có thể thấy rằng:
1Mn0),n1M(h)n(h −≤≤−−−=
và
2
,
2
1M π
±=β
−
=α
(1.18)
Điều này có nghĩa rằng đáp ứng xung h(n) là phản đối xứng (antisymmetric). Chỉ số
đối xứng vẫn là
2
1M −
=α

1M −
=α
không phải là một số nguyên. Đáp
ứng xung được mô tả trong hình 1.8.
1.4.2. Đáp ứng tần số H(e
jω
)
Như vậy, khi tổ hợp hai loại đối xứng và phản đối xứng với M chẵn và M lẻ, ta có
bốn kiểu lọc FIR pha tuyến tính. Đáp ứng tần số của mỗi kiểu có biểu thức và hình
dạng riêng. Để nghiên cứu các đáp ứng pha của các kiểu này, ta viết biểu thức của
H(e
jω
) như sau:
2
1M
,
2
;e)e(H)e(H
)(jj
r
j
−
=α
π
±=β=
αω−βωω
(1.19)
trong đó H
r
(e

=
ω






ω=
∑
(1.20)
trong đó:

( )






−
=
2
1M
h0a
với mẫu ở chính giữa
(1.21)

( )


2
1M −
=α

không phải là một biến nguyên, thì ta có thể chứng tỏ rằng:
( )
( ) 2/1Mj
2/M
1n
j
e
2
1
ncosnb)e(H
−ω−
=
ω

























−ω=ω
2/M
1n
r
2
1
ncosnb)(H
(1.24)
Lưu ý: Tại
π=ω
, ta có
( )
0
2
1
ncosnb)(H
2/M
1n
r

( ) ( )
n1Mhnh −−−=
,
1Mn0 −≤≤
, và
0
2
1M
h =






−
thì ta có thể chứng tỏ:
( )
( )






ω










−
−
= n
2
1M
h2nc
với
2
M
, ,2,1n =
(1.26)
So sánh (1.24) và (1.18), ta có:
( ) ( )
( )
∑
−
=
ω=ω
2/1M
0n
r
nsinncH
(1.27)
Lưu ý: Tại
0=ω

,
( ) ( )
n1Mhnh −−−=
,
1Mn0 −≤≤
, nhưng
2
1M −
=α
không phải là một biến nguyên, thì ta có thể chứng tỏ rằng:
( )






ω






−
−
π
=
ω


(1.28)
trong đó:
( )






−= n
2
M
h2nd
với
2
M
, ,2,1n =
(1.29)
So sánh (1.27) và (1.18), ta có:
( )
∑
=









FIR Type 1
FIR Type 2
FIR Type 3
FIR Type 4
THIẾT KẾ BỘ LỌC FIR BẰNG PHƯƠNG PHÁP LẤY MẪU TẦN SỐ
(FREQUENCY SAMPLING DESIGN TECHNIQUES)
11
DSP – Thiết kế bộ lọc FIR bằng phương pháp lấy mẫu tần số
Cách tiếp cập này đưa trên sự kiện là hàm truyền hệ thống H(z) có thể thu được từ
các mẫu H(k) của đáp ứng tần số H(e
jw
).
Kỹ thuật thiết kế này rất phù hợp với cấu trúc lấy mẫu tần số đã khảo sát trong
chương 5.

Đáp ứng Pha đối với Type 1 & 2
Đáp ứng Pha đối với Type 3 & 4
Đáp ứng thực tế là nội suy của các mẫu được cho bởi12
DSP – Thiết kế bộ lọc FIR bằng phương pháp lấy mẫu tần số
PHẦN 2. THIẾT KẾ BỘ LỌC
1. BỘ LỌC THÔNG THẤP – LÊ PHƯỚC SINH
1.1.Bộ lọc thông thấp lý tưởng:
Bộ lọc số lý tưởng có đặc tính biên dộ tần số dạng chữ nhật:
Trên thực tế không thể xây dựng được bộ lọc số có đặc tính biên độ tần số như vậy,
tuy nhiên các bộ lọc số lý tưởng là cơ sở để phân tích và tổng hợp các bộ lọc số thực
tế.
Bộ lọc thông thấp lý tưởng có đặc tính biên độ tần số khi như sau:

M=20;
alpha=(M-1)/2;
l=0:M-1;
wl=(2*pi/M)*l;
Hrs=[1,1,1,zeros(1,15),1,1];
Hdr=[1,1,0,0];
wdl=[0,0.25,0.3,2];
k1=0:floor((M-1)/2);
k2=(floor((M-1)/2)+1):M-1;
angH=[-alpha*(2*pi)/M*k1,alpha*(2*pi)/M*(M-k2)];
H=Hrs.*exp(j*angH);
h=real(ifft(H,M));
[db,mag,pha,grd,w]=freqz_m(h,1);
[Hr,w,a,L]=Hr_Type2(h);
%subplot(1,1,1)
subplot(2,2,1); plot(wl/pi, Hrs,'o',wdl,Hdr);
axis([0,1,-0.1,1.1]); title('Frequency Samples: M=20')
14
DSP – Thiết kế bộ lọc FIR bằng phương pháp lấy mẫu tần số
xlabel('frequency in pi units'); ylabel('Hr(k)')
subplot(2,2,2);
stem(l,h);axis([-1,M,-0.1,0.3])
title('Impulse Response');
xlabel('n');
ylabel('h(n)');
subplot(2,2,3);
plot(w/pi,Hr,wl/pi,Hrs,'o')
axis([0,1,-0.2,1.2]); title('Amplitude Response')
xlabel('frequency in pi units'); ylabel('Hr(w)')
subplot(2,2,4); plot(w/pi,db); axis([0,1,-60,10]);grid

function [db,mag,pha,grd,w]=freqz_m(b,a);
% Modified version of freqz subroutine
% [db,mag,pha,grd,w]=freqz_m(b,a);
% db=Relative magnitude in dB computed over 0 to pi radians
% mag=absolute magnitude computed over 0 to pi radians
% grd= Group delay over 0 to pi radians
% w=501 frequency samples between 0 to pi radians
% b=numerator polynomial of H(z) (for FIR: a=h)
% a=demonitor polynomial of H(z) (for FIR: a=[1])
[H,w]=freqz(b,a,1000,'whole');
16
DSP – Thiết kế bộ lọc FIR bằng phương pháp lấy mẫu tần số
H=(H(1:1:501))';w=(w(1:1:501))';
mag=abs(H);
db=20*log10((mag+eps)/max(mag));
pha=angle(H);
grd=grpdelay(b,a,w);
Chương trình thực hiện.
M=33; alpha=(M-1)/2;
l=0:M-1;
w1=(2*pi/M)*l;
T1=0.12; T2=0.8;
Hrs=[zeros(1,11),T1,T2,ones(1,8),T2,T1,zeros(1,10)];
Hdr=[0,0,1,1];
wdl=[0,0.6,0.8,1];
k1=0:floor((M-1)/2);
k2=floor((M-1)/2+1):(M-1);
angH=[-alpha*(2*pi)/M*k1,alpha*(2*pi)/M*(M-k2)];
H=Hrs.*exp(j*angH);
h=real(ifft(H,M));

dBR
dBA
Ss
pp
pp
SS
60;8,0
1;65,0
1;35,0
60;2,0
2
2
1
1
==
==
==
==
πω
πω
πω
πω
- Với A
S
=60dB chọn cửa sổ Blackman với hàm cửa sổ như sau:






ππ
- Tính M với :
πππωωω
πππωωω
15,065.08,0
15,02,035,0
222
111
=−=−=∆
=−=−=∆
PS
SP
πωω
15.0),(
21
=∆∆=→ MINTW
→=+=+=⇒ 3.741
15,0
11
1
11
π
ππ
TW
M
chọn M=75
- Tìm tần số cắt:
( )
( ) ( )
πππωωω

2
22
1
2
1
21
2
ω
ω
π
ω
ω
ω
π
ω
π
ωω
ω
ω
ω
ω
ω
ω
−=







p1
ω
đó là tại k=7:
7
40
2
2.0
1
π
πω
==
p
19
DSP – Thiết kế bộ lọc FIR bằng phương pháp lấy mẫu tần số
Và mẫu tiếp theo tại
2s
ω
đó là tại k=13:
13
40
2
2.0
2
π
πω
==
p
Và mẫu tiếp theo tại
p2
ω

5.19
2
140
=
−
=
α
với giá trị này thì thuộc bộ đáp ứng pha tuyến
Tính Type-2, từ [7.37] chúng ta có:





≤≤−+
≤≤−=−
=∠
3920),40(95.1
190,5.19
20
2
5.19
)(
kk
kkk
kH
π
π
π
Tìm đáp ứng xung h(n)

subplot(1,1,1)
subplot(2,2,1);plot(wl(1:21)/pi,Hrs(1:21),'o',wdl,Hdr);
axis([0,1,-0.1,1.1]); title('Frequency Samples : M=40')
xlabel(' '); ylabel('Hr(k)')
set(gca,'XTickMode','manual','XTick',[0,0.2,0.35,0.65,0.8,1])
set(gca,'YTickMode','manual','YTick',[0,1]); grid
subplot(2,2,2); stem(l,h); axis([-1,M,-0.4,0.4]);grid
title('Impulse Response');ylabel('h(n)');text(M+1,-0.1,'n')
subplot(2,2,3); plot(ww/pi,Hr,wl(1:21)/pi,Hrs(1:21),'o');
axis([0,1,-0.2,1.2]); title('DAP UNG BIEN DO')
xlabel('Frequency in pi units'); ylabel('Hr(w)')
set(gca,'XTickMode','manual','XTick',[0,0.2,0.35,0.65,0.8,11])
set(gca,'YTickMode','manual','YTick',[0,1]); grid
subplot(2,2,4);plot(w/pi,db); axis([0,1,-100,10]); grid
title('Magnitude Response '); xlabel('Frequency in pi units');
ylabel('dB');
set(gca,'XTickMode','Manual','XTick',[0,0.2,0.35,0.65,0.8,1]);
set(gca,'YTickMode','Manual','YTick',[-16;0]);
set(gca,'YTickLabelMode','manual','YTickLabels',['16';' 0'])
* Phương pháp thiết kế tối ưu:
% Phuong phap toi uu
% ws1=0.2pi, wp1=0.35pi, wp2=0.65pi, ws2=0.8pi, Rp=1dB,
As=60dB
% T2 = 0.59417456, T1=0.109021
M = 40; alpha = (M-1)/2; l = 0:M-1; wl = (2*pi/M)*l;
T1 = 0.109021; T2 = 0.59417456;
Hrs =
[zeros(1,5),T1,T2,ones(1,7),T2,T1,zeros(1,9),T1,T2,ones(1,7),T2,T1,zeros
(1,4)];
Hdr = [0,0,1,1,0,0]; wdl = [0,0.2,0.35,0.65,0.8,1];

Hình 4.3: Mô tả đáp ứng sau khi thiết kế
Nhận xét:
- Lỗi xấp xỉ là hiệu của đáp ứng lý tưởng và đáp ứng thực tế bằng không tại các
tần số được lấy mẫu.
22
DSP – Thiết kế bộ lọc FIR bằng phương pháp lấy mẫu tần số
- Lỗi xấp xỉ ở tất cả các tần số khác nhau phụ thuộc vào hình dạng của đáp ứng
tần số lý tưởng; nghĩa là, đáp ứng tần số lý tưởng càng sắc nét thì lỗi xấp xỉ càng
lớn.
- Lỗi càng lớn khi ở gần cạnh dải và càng bé khi ở bên trong dải.
Quan sát độ thị hình 4.3 thì độ suy hao dải chắn tối thiểu khoảng 16dB không
thỏa mãn với giá trị As=60dB đã cho. Do vậy phương pháp xứ lý này không chấp
nhận được.
Để có được hệ số suy giảm tốt hơn, chúng ta phải tăng M tạo ra các mẫu tự do ở
dải chuyển tiếp nghĩa là, chúng ta thay đổi các giá trị của chúng để thu được hệ số suy
giảm lớn nhất đối với M và độ rộng dải chuyển tiếp đã cho. Trong thực tế thì trong dải
chuyến tiếp người ta chọn một hoặc hai mẫu. Đây được gọi là phương pháp thiết kế tối
ưu.
* Phương pháp thiết kế tối ưu:
Hình 4.4: Mô tả đáp ứng sau khi thiết kế
Nhận xét:
Độ suy hao dải chẵn tối thiểu bây giờ là 60dB. So với giá thiết thì có thể chấp
nhận được. Do đó chấp nhận bộ lọc FIR pha tuyến tính mà đáp ứng tấn số của nó độ
gợn song của biên độ là hằng số.
23
DSP – Thiết kế bộ lọc FIR bằng phương pháp lấy mẫu tần số
4. THIẾT KẾ BỘ LỌC FIR CHẮN DẢI - VŨ BẢO TOÀN
BÀI TOÁN
4.1. Thiết kế bằng phương pháp lấy mẫu tần số:
Các bước tính toán thiết kế lọc FIR chắn dải với các thông số sau:

222
111
=−=−=∆
=−=−=∆
sp
ps
Tần số cắt:
( )
( ) ( )
πππωωω
πππωωω
65,02/6,07,02/
35,02/)4,03,0(2/
222
111
=+=+=
=+=+=
SPC
SPC
Chúng ta chọn M=63 do đó chúng ta có mẫu tần số tại
p1
ω
đó là tại k = 9:
9
63
2
3.0
1
π
πω

Và mẫu tiếp theo tại
p2
ω
đó là tại k= 22:
22
63
2
7.0
2
π
πω
==
p
Do đó chúng ta có 9 mẫu trong stopband [
p1
0
ωω
≤≤
] và 6 mẫu trong
passband [
ss 21
ωωω
≤≤
] và 9 mẫu trong stopband [
πωω
≤≤
p2
]. Chúng ta có:
( )
KH

π
π
Tìm đáp ứng xung h(n)
4.3. Chương trình Matlab
Các hàm sử dụng
Hàm Hr_Type1(h):Tính đáp biên độ của bộ lọc tần số loại 1
function [Hr,w,a,L] = Hr_Type1(h);
% Computes Amplitude response Hr(w) of a Type-1 LP FIR filter
% [Hr,w,a,L] = Hr_Type1(h)
% Hr = Amplitude Response
% w = 500 frequencies between [0 pi] over which Hr is computed
% a = Type-1 LP filter coefficients
% L = Order of Hr
% h = Type-1 LP filter impulse response
M = length(h);
L = (M-1)/2;
a = [h(L+1) 2*h(L:-1:1)]; % 1x(L+1) row vector
n = [0:1:L]; % (L+1)x1 column vector
w = [0:1:500]'*pi/500;
Hr = cos(w*n)*a';
freqz_m(b,a): Tính đáp ứng tần số của bộ lọc FIR có đáp ứng xung h.
function [db,mag,pha,grd,w]=freqz_m(b,a);
% Modified version of freqz subroutine
%
25