Filtry cyfrowe
Procesory sygnałowe (DSP),
układy programowalne
x(n) y(n)
Filtr cyfrowy
h(n)
h(n)
 odpowiedz
impulsowa
x(n) y(n)
1
y(n) = x(n) " h(n)
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Filtry cyfrowe
Po co filtrujemy sygnały?
Aby uzyskać:
redukcję zakłóceń sygnału
(np. zakłóceń od sieci energetycznej)
(np. zakłóceń od sieci energetycznej)
zmianę charakterystyki widmowej sygnału
(preemfaza, deemfaza)
wyodrębnienie zadanych składowych sygnału
spośród jego innych składowych (detekcja)
2
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Charakterystyki
częstotliwościowe - filtry idealne
A(f)
A(f)
górnoprzepustowy
dolnoprzepustowy
f
f
f
A(f)
A(f)
pasmowoprzepustowy
pasmowozaporowy
f
f
3
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Charakterystyki
częstotliwościowe filtrów
pasmo
A(f)
przejściowe
Filtr dolno-przepustowy
pasmo
(np. filtr anty-alisingowy,
zaporowe
redukcja zakłóceń)
pasmo
przepustowe
f
fp fs
A(f)
Filtr górno-przepustowy
(np. preemfaza)
4
f
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Charakterystyki
częstotliwościowe filtrów
A(f)
Filtr środkowo-przepustowy
(np. detekcja cech sygnału)
f
f
A(f)
Filtr środkowo-zaporowy
(np. redukcja zakłóceń od sieci
energetycznej)
50 Hz
f
5
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Charakterystyki
częstotliwościowe filtrów - przykłady
A(f)
Filtr dolnoprzepustowy
(redukcja zakłóceń)
f
1500 1500
1000 1000
500 500
0 V 0 V
0 0
6
200 400 600 800 1000 200 400 600 800 1000
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Charakterystyki
częstotliwościowe filtrów - przykłady
A(f)
Filtr górno-przepustowy
(np. usuwanie wartości średniej)
f
100
1500
50
1000
0
-50
500
-100
0 V
0 V
0
-150
0 200 400 600 800 1000
200 400 600 800 1000
7
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Charakterystyki
częstotliwościowe filtrów - przykłady
A(f)
Filtr środkowo-przepustowy
(np. detekcja cech sygnału)
f
f
6000
1500
4000
1000
2000
500
0
0 V
0
-2000
200 400 600 800 1000
0 200 400 600 800 1000
8
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Charakterystyki
częstotliwościowe filtrów - przykłady
A(f)
Filtr środkowo-zaporowy
(np. redukcja zakłóceń o
zadanej częstotliwości)
50 Hz
50 Hz
f
f
1600
1400
1500
1300
1400
1200
1300
1200
1100
1100
1000
1000
900
900
0 200 400 600 800 1000
0 200 400 600 800 1000
9
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Zastosowania filtrów cyfrowych w
przetwarzaniu elektrokardiogramu
donoprzepustowe (redukcja zakłóceń o
częstotliwościach radiowych, aktywności mięśni
szkieletowych)
górnoprzepustowe (eliminacja pełzania linii
górnoprzepustowe (eliminacja pełzania linii
izoelektrycznej, fg=0.5 Hz, zob.  ecg_mit.mat )
pasmowoprzepustowe (wydzielanie
składowych sygnału EKG, np. fali P, T, QRS)
pasmowozaporowe (redukcja zakłóceń od sieci
energetycznej, f=50 Hz)
10
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Rodzaje filtrów
A
B
Filtr
Filtr
?
?
C
- Dolnoprzepustowy
- Pasmowo-przepustowy
- Górnoprzepustowy
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Filtry cyfrowe  SOI i NOI
1
Filtry dzielimy również na:
0.8
0.6
filtry o skończonej
0.4
odpowiedzi impulsowej (SOI)
0.2
0
tzw. filtry nierekursywne
tzw. filtry nierekursywne
0 5 10 15 20
1
0.8
filtry o nieskończonej
0.6
odpowiedzi impulsowej (NOI)
0.4
0.2
tzw. filtry rekursywne
0
0 5 10 15 20
12
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Filtry cyfrowe  SOI i NOI
Elementy wykonawcze potrzebne do realizacji
filtru cyfrowego:
Mnożenie Opóz
nienie
Sumator
Sumator
przez stałą jednostkowe
u1(k)
y(k)
y(k)
c
u(k)
u(k)
z-1 y(k)
u2(k)
y(k)= cu(k) y(k)= u(k - 1)
y(k)= u1(k)+ u2(k)
13
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Idea filtracji cyfrowej
współczynniki
współczynniki
ruchomej średniej
autoregresji
x(k) y(k)
a1=1
b
0
z-1 z-1
a
a
b
b
1 2
1 2
x(k-1) y(k-1)
z-1 z-1
b a
2 3
x(k-2) y(k-2)
z-1 z-1
bM aN
14
x(k-M) y(k-N)
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Filtr o Skończonej odpowiedzi
Impulsowej(SOI) - przykład
LP:
b0 =0.5
x(k)
y(k)= 0.5x(k)+ 0.5x(k -1)
z
z-1
"t
"t
b1 =0.5
x(k-1)
HP:
y(k)= 0.5x(k)- 0.5x(k -1)
BP: ?
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Filtr o Nieskończonej odpowiedzi
Impulsowej (NOI) - przykład
y(k)= a1y(k -1)+ x(k)
x(k) y(k)
a1 "t
a1y(k-1)
y(k-1)
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Równanie różnicowe filtru
M N
y(k)=
"b(n)x(k - n)-"a(n)y(k - n)
n=0 n=1
Z uwagi na implementację programową filtrów w
Z uwagi na implementację programową filtrów w
programie Matlab (zob. funkcję  filter ) stosujemy
następujące indeksowanie współczynników filtru:
M N
*
1=
a(1)y(k)=
"b(n)x(k - n +1)-"a(n)y(k - n +1)
n=1 n=2
17
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Równanie różnicowe filtru
M N
y(k)=
*
"b(n)x(k - n +1)-"a(n)y(k - n +1)
n=1 n=2
Jeżeli wszystkie współczynniki a(n) są zerowe to
Jeżeli wszystkie współczynniki a(n) są zerowe to
równanie różnicowe opisuje filtr cyfrowy SOI, w
przeciwnym przypadku filtr NOI
SOI  ang. Finite Impulse Response (FIR)
NOI  ang. Infinite Impulse Response (IIR)
18
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Równanie różnicowe filtru
Zwróćmy uwagę na indeksowanie
współczynników filtru i próbek sygnału:
y(k)= b(1)x(k) + b(2)x(k -1)+ b(3)x(k - 2)+K
y(k)= b(1)x(k) + b(2)x(k -1)+ b(3)x(k - 2)+K
*
*
- a(2)y(k -1)- a(3)y(k - 2)- a(4)y(k - 3)-K
Jest to równanie filtru przyczynowego, & ...
19
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Prosty filtr górnoprzepustowy
(freqz)
0
-20
-40
h=[0.5 -0.5]
-60
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
Normalized Frequency ( rad/sample)
�Ą
�Ą
100
50
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
Normalized Frequency ( rad/sample)
�Ą
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Magnitude (dB)
Phase (degrees)
Prosty filtr dolnoprzepustowy
(freqz)
-3dB
0
-20
-40
h=[0.5 -0.5]
-60
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
Normalized Frequency ( rad/sample)
�Ą
�Ą
0
-50
-100
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
Normalized Frequency ( rad/sample)
�Ą
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Magnitude (dB)
Phase (degrees)
Prosty przykład filtru SOI
Filtr o ruchomej średniej:
k
1 1
y(k)= [x(k - 4)+ x(k - 3)+ x(k - 2)+ x(k -1)+ x(k)]=
"x(n)
5 5
n=k -4
Określ wektory współczynników a i b dla tego filtru (wzór
*).
Odpowiedz
:
a=[1];
b=[0.2 0.2 0.2 0.2 0.2];
Zatem jest to filtr SOI   kiepski filtr dolnoprzepustowy.
22
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Prosty przykład filtru SOI
Jaka jest odpowiedz
impulsowa tego filtru?
k
1 1
y(k)= [x(k - 4)+ x(k - 3)+ x(k - 2)+ x(k -1)+ x(k)]=
"x(n)
5 5
n=k -4
Odpowiedz
:
h=[0.2 0.2 0.2 0.2 0.2];
Wniosek:
Wniosek: dla filtru SOI współczynniki filtru i jego
!
!
odpowiedz
impulsowa to jest to samo!
23
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Prosty przykład filtru SOI
Ćwiczenie:
Ćwiczenie:
Korzystając z instrukcji:
- conv
- filter
- filter
Wyznacz odpowiedz
zdefiniowanego filtru SOI na
pobudzenie sygnałem:
x=[10 50 24 60 37 77 89 22 63 9 52 31 49 54 28 14];
Na tym samym wykresie wyświetl sygnał oraz
otrzymane sygnały po filtracji. Porównaj długość
24
sygnałów wynikowych?
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Prosty przykład filtru SOI
Ćwiczenie:
Ćwiczenie:
Wyznacz charakterystykę amplitudową i fazową
zadanego filtru.
Podpowiedz
:
- utwórz wektor 512 elementowy uzupełniając wektor
odpowiedzi impulsowej 507 zerami,
- oblicz transformacje Fouriera tego wektora (fft)
- wyznacz charakterystyki amplitudowe (abs) i fazowe
(angle) filtru przyjmując, że zastosowano fs=100 Hz.
25
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Prosty przykład filtru SOI
Charakterystyka fazowa Charakterystyka amplitudowa
Charakterystyka fazowa
Charakterystyka amplitudowa
3
1
2
0.8
1
tzw. zero filtru
0.6
0
0.4
-1
-1
0.2
-2
-3 0
0 20 40 60 80 100 0 20 40 60 80 100
f [Hz] f [Hz]
0.2
0.15 %MATLAB
Odpowiedz

freqz(b,a)
0.1
impulsowa filtru
freqz(b,a,N,fs)
0.05
26
0
0 5 10 15 20
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
radiany
Amplituda
Prosty przykład filtru SOI
Ćwiczenie cd:
Ćwiczenie cd:
Zbadaj amplitudy i przesunięcia fazowe przebiegów
wyjściowych filtru gdy na jego wejście podano:
przebieg sinusoidalny f=1 Hz
przebieg sinusoidalny f=1 Hz
przebieg sinusoidalny f=5 Hz
przebieg sinusoidalny f=20 Hz
Potwierdz
, dla zadanego sygnału wejściowego filtru,
że filtracja (splot) w dziedzinie czasu odpowiada
mnożeniu widma odpowiedzi impulsowej filtru i widma
27
sygnału.
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Prosty przykład filtru SOI
1 Hz
1
0.8
0.6
0.4
0.2
0
-0.2
-0.4
-0.6
-0.8
-1
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Prosty przykład filtru SOI
5 Hz
1
0.8
0.6
0.4
0.2
0
-0.2
-0.4
-0.6
-0.8
-1
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Prosty przykład filtru SOI
20 Hz
1
0.8
0.6
0.4
0.2
0
-0.2
-0.4
-0.6
-0.8
-1
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Prosty przykład filtru SOI
input signal
2
1.5
1
0.5
0
-0.5
-1
-1.5
-2
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Prosty przykład filtru SOI
xwy = a sin(10Ąt +ą)+ b sin(60Ąt + � )
xwe = sin(10Ąt)+ sin(60Ąt)
fs = 200Hz
a,b,ą ,� = ?
2
1.5
1.5
1
1
1
0.5
0.5
0
0
-0.5
-0.5
-1
-1
-1.5
-2
-1.5
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
[sek]
[sek]
tzw. liniowe zniekształcenia amplitudowe
32
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Wej
Wyj
Przekształcenie z
Definicja:
Definicja: j�t0
x(t - t0)"! X (�)e
Przekształceniem z ciągu próbek:
j�
{x(0), x(T), x(2T),K, x(nT),K}
z = e
j�(n-1) j�n
j�(n-1) j�n
jest wyrażenie:
jest wyrażenie:
e = e- j�e
e = e- j�e
N
X (z)= x(0)+ x(T )z-1 + x(2T)z-2 +K+ x(NT)z-N = x(nT)z-n
"
n=0
Przekształcenie z pełni ważną rolę w cyfrowym
przetwarzaniu sygnałów, np. dla impulsu
jednostkowego X(z)=1
33
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Przekształcenie z - analogia do DFT
Przekształcenie z:
" " "
-n
-n j� -n jn�
X (z)= (re ) =
"x(n)z = "x(n) "x(n)r e-
-" -" -"
Zatem dla r=1 równanie powyższe,
|z|=1
tj. przekształcenie z jest równoważne
r=1
szeregowi Fouriera.
Jeżeli x(n) reprezentuje odpowiedz

impulsową filtru, to X(z) dla |z|=1 jest
charakterystyką częstotliwościową filtru.
z=ej�
34
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Przekształcenie z - analogia do DFT
Przekształcenie z:
N -1
X(z)= x(n)z-n
"
n=0
2Ą
�ł �ł
jk
jk
�ł �ł
�ł �ł
j�
j�
N
�ł łł
Przekształcenie dyskretne Fouriera:
Przekształcenie dyskretne Fouriera:
z = e = e
= =
2Ą 2Ą
�ł �ł
jk�ł �ł - j k�ł �ł n
�ł �ł �ł �ł
N -1
�ł �ł
j� N N
�ł łł �ł łł
X(e )= X = X (k)= x(n)e
"
�łe �ł
�ł �ł n=0
�ł łł
indeks czasu
�=2Ąf/fs
indeks częstotliwości
35
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Przekształcenie z
Równanie różnicowe filtru SOI:
y(nT)= x(nT)+ 2x(nT -T )+ 3x(nT - 2T)
Można wyrazić w dziedzinie przekształcenia z w postaci:
Y(z)= X (z)+ 2X (z)z-1 + 3X (z)z-2
i otrzymać funkcję przejścia filtru:
Y(z)
H(z)= = 1+ 2z-1 + 3z-2
X (z)
36
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Przekształcenie z
Podobnie, ogólne równanie różnicowe filtru cyfrowego:
M N
y(k)=
"b(n)x(k - n)-"a(n)y(k - n)
n=0 n=1
w dziedzinie przekształcenia z można zapisać w postaci:
M
zera filtru
"b z-n
n
Y(z)
n=0
H(z)= =
N
X (z)
bieguny filtru
1+
"a z-n
n
n=1
37
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Płaszczyzna z
Im(z)
Zmienną z definiuje się:
z=j
radiany
j�
z = e
Ą
na okres
2
r=1
r=1
f z=1
f z=1
z=-1
z=-1
� =Ą
� =Ą
� = 2Ą
� = 2Ą
0
0
fs
�=2Ą
fs
f = fs Re(z)
f =
2
3
pulsacja unormowana
2Ą
względem fs
z=-j
Filtr jest stabilny gdy bieguny filtru leżą wewnątrz
38
okręgu jednostkowego.
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Płaszczyzna z
%MATLAB
zplane(0.2*ones(1,5),1)
0.4Ą
0.8Ą
Charakterystyka amplitudowa
Charakterystyka amplitudowa
1
1
0.8
0.5
tzw. zero filtru
0.6
0
0.4
-0.5
0.2
-1
0
-1 -0.5 0 0.5 1
0 20 40 60 80 100
f [Hz]
Real part
39
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Amplituda
Imaginary part
Podstawowe filtry SOI
Ćwiczenie:
Ćwiczenie:
Wyznacz odpowiedz
impulsową i analitycznie wyznacz
charakterystykę amplitudową filtrów:
H(z)= 0.5(1+ z-1)
~ pierwsza pochodna sygnału
H(z)=1- z-1
Podpowiedz
: przyjmij T=1 i zastosuj podstawienie:
j�T
z = e
40
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Podstawowe filtry SOI
Ćwiczenie rozwiązanie :
Ćwiczenie rozwiązanie :
j�
H(z)= 0.5(1+ z-1)= 0.5(1+ e- j�)= H(e )
j� j� 2
H(e ) = 0.5(1+ e- j�) = 0.5e- j� 2(e + e- j� 2) =
� �
charakterystyka
0.5 e- j� 2 2cos = cos
amplitudowa
2 2
� �
�ł
j�
charakterystyka
arg H(e )= arg�łe- j� 2 cos = -
�ł �ł
fazowa
2 2
41
�ł łł
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Podstawowe filtry SOI
0
-20
-40
-60
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
Normalized frequency (Nyquist == 1)
0
%MATLAB
-20
freqz(0.5*[1 1],1);
-40
-60
-80
-100
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
Normalized frequency (Nyquist == 1)
42
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Magnitude Response (dB)
Phase (degrees)
Warunek liniowej fazy dla filtrów SOI
Można projektować filtry SOI tak by zachować
liniową charakterystykę fazową. Zapewniają to
warunki symetrii dla współczynników filtru:
h(M - k)= h(k)
h(M - k)= h(k)
dla
k = 0,1,K, M
dla M parzystych lub nieparzystych
h(M - k)= -h(k)
Jaka jest korzyść z liniowej fazy filtru?
Opóz
nienie sygnału na wyjściu filtru jest niezależne od
!
!
jego częstotliwości (zob. pierwsze ćwiczenie nt. filtrów).
43
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Filtr o liniowej fazie
Ą/4  opóz
nienie fazowe
s(t)=sin(�t)
td  opóz
nienie czasowe
s(t)=sin(2�t)
Ą/2  opóz
nienie czasowe
Filtry o liniowej fazie opóz
niają wszystkie
częstotliwości o ten sam czas!
44
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Warunek liniowej fazy dla filtrów SOI
h(M - k)= h(k)
dla
k = 0,1,K, M
h(M - k)= -h(k)
czas dyskretny
czas dyskretny
czas dyskretny
czas dyskretny
45
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Zniekształcenia fazowe
Ą Ą
�ł �ł
xwy = sin�ł10Ąt + + sin�ł20Ąt +
�ł �ł �ł �ł
xwe = sin(10Ąt)+ sin(20Ąt)
2 2
�ł łł �ł łł
2 2
1.5
1.5
1
1
0.5
0.5
0
0
-0.5
-0.5
-1
-1
-1.5
-2 -1.5
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
[sek] [sek]
tzw. nieliniowe zniekształcenia fazowe
(np. duże zniekształcenia sygnałów akustycznych)
46
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Wyj
Wyj
Zniekształcenia fazowe - przykład
30
20
20
10
0
-10
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
Normalized frequency (Nyquist == 1)
0
-20
%Matlab -40
-60
x=wavread( glos_meski.wav );
-80
y=filter(1,[1 -0.95],x);
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
Normalized frequency (Nyquist == 1)
sound(y,fs)
47
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Magnitude Response (dB)
Magnitude Response (dB)
Phase (degrees)
Jak przekształcić dolnoprzepustowy
filtr SOI w filtr górno przepustowy?
Ćwiczenie:
Ćwiczenie:
Zastosuj następujące przekształcenie dla współczynników
bH przykładowego dolnoprzepustowego filtru SOI.
n
bH = (-1) bL(n)
Zbadaj charakterystykę otrzymanego filtru i porównaj ją z
prototypowym filtrem dolnoprzepustowym.
Zobacz tw. o modulacji dla przekształcenia Fouriera
48
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Inne przykłady filtrów SOI
Ćwiczenie cd:
Ćwiczenie cd:
Zbadaj charakterystyki filtrów o odpowiedzi impulsowej:
h=0.25*[1 2 1] , tzw. filtr von Hanna
h=0.25*[1 2 1] , tzw. filtr von Hanna
h=[1 -2 1] , ~ druga pochodna sygnału
49
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Projektowanie filtrów SOI
metodą rozmieszczania zer filtru
Równanie różnicowe filtru SOI:
Y(z)= b(0)X(z)+ b(1)X(z)z-1 +K+ b(n)X(z)z-n
Transmitancja operatorowa filtru SOI:
Transmitancja operatorowa filtru SOI:
H(z)= b(0)+ b(1)z-1 +K+ b(n)z-n =
�ł �ł
b(1) b(n)�ł =
= b(0)z-n�ł zn + zn-1 +K+
�ł
b(0) b(0)�ł
�ł łł
= b(0)z-n(z - z1)(z - z2)K(z - zn)
50
zera filtru
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Projektowanie filtrów SOI
metodą rozmieszczania zer filtru
Pierwiastki równania:
H(z)= b(0)z-n(z - z1)(z - z2)K(z - zn)= 0
to liczby z1, z , & zn. Ponieważ współczynniki
to liczby z , z2, & z . Ponieważ współczynniki
równania są rzeczywiste to zera tworzą pary
zespolone sprzężone albo są rzeczywiste.
Charakterystyki częstotliwościowe filtrów SOI
można kształtować odpowiednio rozmieszczając
zera filtru.
51
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Projektowanie filtrów SOI
metodą rozmieszczania zer filtru
Przykład:
Skonstruować filtr SOI, który  zeruje pulsację Ą/4
rad/okres.
1
0.8
jĄ / 4
( )= (0) ( - )( - )
H(z)= b(0)z-2(z - e )(z - e- jĄ / 4)
0.6
0.6
0.4
0.2
Dla b(0)=1:
0
-0.2
jĄ / 4
-0.4
H(z)= z-2(z - e )(z - e- jĄ / 4)=
-0.6
Ą
-0.8
= z-2�ł z2 - 2z cos +1�ł =1- 2z-1 + z-2
�ł �ł
-1
4
�ł łł
-1 -0.5 0 0.5 1
Real part
Równanie różnicowe filtru:
y(k)= x(k)- 2x(k -1)+ x(k - 2)
52
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Imaginary part
Projektowanie filtrów SOI
metodą rozmieszczania zer filtru
Przykład:
kĄ
Sprawdzenie dla sygnału:
x(k)= sin
4
k : 0 1 2 3 4 5 6
k : 0 1 2 3 4 5 6
x(k): 0 1/ 2 1 1/ 2 0 -1/ 2 -1
y(k): 0 1/ 2 0 0 0 0 0
Równanie różnicowe filtru:
y(k)= x(k)- 2x(k -1)+ x(k - 2)
53
Sprawdz
przebieg charakterystyki widmowej: freqz
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Projektowanie filtrów SOI
metodą okien czasowych
Chcemy zaprojektować idealny
A(�)
filtr dolnoprzepustowy.
Otrzymujemy nierealizowalną,
nieskończoną w czasie
nieskończoną w czasie
charakterystykę odpowiedzi
�0 �
impulsowej:
+�0
1
- j�n
h(n)= d� = ?
+"e
2Ą
-�0
Należy ograniczyć czas trwania tej odpowiedzi.
54
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Projektowanie filtrów SOI
metodą okien czasowych
Zastosowanie okna czasowego ograniczającego czas
trwania tej odpowiedzi pozwala uzyskać filtr realizowalny
fizycznie,
np. dla filtru dolnoprzepustowego o częstotliwości odcięcia
0.4Ą rad/s i odpowiedzi impulsowej ograniczonej do 51
próbek:
próbek:
b=0.4*sinc(0.4*(-25:25));
%zobacz równie\
(-100:100)
uzyskuję się charakterystykę:
[H,f]=freqz(b,1,512,2);
plot(f,abs(H)),grid
55
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Projektowanie filtrów SOI
metodą okien czasowych
1.4
1.2
1
0.8
tzw. efekt Gibbsa
0.6
~9% amplitudy impulsu
~9% amplitudy impulsu
0.4
0.2
0
f
0 0.2 0.4 0.6 0.8 1
Efekt Gibbsa można zredukować stosując zamiast
okna prostokątnego wycinającego odpowiedz

impulsową, okno o kształcie podobnym do funkcji
56
Gaussa, np. okno Hamminga
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Projektowanie filtrów SOI
metodą okien czasowych
1
Okno Hamminga
0.8
0.6
%MATLAB
0.4
0.4
b=b.*hamming(51) ;
0.2
[H,f]=freqz(b,1,512,2);
f
0
10 20 30 40 50
plot(f,abs(H)),grid
rząd filtru
W programie Matlab opisaną procedurę projektowania
filtrów implementuje instrukcja syntezy filtru FIR  fir1
57
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Zadanie dla dociekliwych
Wyjaśnić skąd się bierze tzw. efekt Gibbsa
efekt Gibbsa.
Pokazać, że maksymalna amplituda oscylacji
w tym efekcie jest niezależna od czasu
w tym efekcie jest niezależna od czasu
trwania okna prostokątnego wycinającego
sygnał.
58
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Prosty filtr NOI
Rozważmy prosty filtr NOI: y(n)= �y(n -1)+ x(n)
Y(z) 1
H(z)= =
X (z) 1- �z
X (z) 1- �z-1
zero z=0
1 z z
H(z)= =
1- �z-1 z z - �
biegun z=�
a(1) a(2)=-�
59
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Prosty filtr NOI
10
1
5
0.5
0
� = 0.5
� = 0.5
-5
-5
0
0 0.2 0.4 0.6 0.8 1
Normalized frequency (Nyquist == 1)
0
-0.5
-10
-1
-20
-1 -0.5 0 0.5 1
Real part
-30
0 0.2 0.4 0.6 0.8 1
Normalized frequency (Nyquist == 1)
60
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Magnitude Response (dB)
Magnitude Response (dB)
Imaginary part
Imaginary part
Phase (degrees)
Prosty filtr NOI
5
x 10
7
6
5
1
4
0.8
3
3
�=1.5>1
0.6 �=1.5>1
�=0.5<1
�=0.5<1
2
0.4
1
0.2
0
0
0 10 20 30
0 5 10
61
pł. z
pł. z
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Projektowanie filtrów NOI
Metoda bezpośrednia - aproksymacyjna:
% MATLAB
% [b,a]=yulewalk(n,f,m)
% n  rząd filtru
% f  próbki char. częstotl. z zakresu <0,1>
% f  próbki char. częstotl. z zakresu <0,1>
% m  dyskretne częstotl. z zakresu <0,1>
f = [0 0.6 0.6 1];
m = [1 1 0 0];
[b,a] = yulewalk(8,f,m);
[h,w] = freqz(b,a,128);
plot(f,m,w/pi,abs(h),'--')
62
Nieliniowa faza! Zobacz też  zplane(b,a)
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Projektowanie filtrów NOI
Metoda niezmienności odpowiedzi impulsowej:
Wyznacz odpowiedzi impulsowe
% MATLAB
tych filtrów
%dolnoprzepustowy Butterwotha
[b,a]=butter(5,0.4)
[b,a]=butter(5,0.4)
%pasmowoprzepustowy Czebyszewa typu I
[b,a]=cheby1(4,1,[.4 .7])
%górnoprzepustowy Czebyszewa typu II
[b,a]=cheby2(6,60,.8, high )
%pasmowozaporowy eliptyczny
[b,a] = ellip(3,1,60,[.4 .7], stop );
63
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Porównanie filtrów SOI i NOI
NOI
SOI
z definicji stabilne mogą być niestabilne
łatwe projektowanie bardziej złożone projektowanie
nieliniowa faza
łatwo zapewnić liniową fazę
łatwo zapewnić liniową fazę nieliniowa faza
uzyskanie stromej możliwość uzyskiwania bardzo
charakterystyki jest możliwe stromej charakterystyki
dla dużego rzędu filtru
problemy implementacyjne
z uwagi na skończoną
skończoną dokładność
reprezentacji współczynników dokładność reprezentacji
filtru nie jest dokuczliwa współczynników filtru
64
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Narzędzie projektowania filtrów  sptool
65
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Idea filtru adaptacyjnego
kolejne próbki sygnału
x0= 1
w1
w0
x(n-P+1)
y(n)
w2
...
"
"
x(n-1)
wP
�(n)
x(n)
"
d(n)
Funkcja celu:
E[d(n)  wTx] = [(d(n) - y(n))2] = E[�2] min
66
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Adaptacyjna redukcja zakłóceń
z
ródło
sygnału
^
e(t) = x(t)
s(t) = x(t) + v(t)
+
Ł
_
^
n(t) v(t)
Filtr
adaptacyjny
z
ródło
(tzw. wejście
zakłóceń
Reguła
odniesienia)
adaptacji wag
"w(t)= -�"e(t)
67
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Adaptacyjna redukcja zakłóceń
Minimum wielkości e(t) odpowiada najskuteczniejszej
(zgodnie z minimum błędu średniokwadratowego)
redukcji zakłóceń:
e(t) = s(t) - v(t) = x(t) + v(t) - v(t)
e(t) = s(t) - v(t) = x(t) + v(t) - v(t)
Ć Ć
Ć Ć
E[e2] = E[x2]+2E[x(v -v)]+ E[(v -v)2] H" E[x2]
Ć Ć
c.n.d.
"
= "
sygnał i zakłócenie
są nieskorelowane
68
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Zastosowania filtracji adaptacyjnej
Filtry adaptacyjne są stosowane głównie do
filtracji sygnałów niestacjonarnych, np.:
w adaptacyjnej redukcji zakłóceń mierzonego sygnału
od sieci energetycznej oraz redukcji zakłóceń od
od sieci energetycznej oraz redukcji zakłóceń od
elektronarzędzi chirurgicznych (f~120 Hz)
do redukcji energii sygnału EKG matki przy pomiarze
EKG płodu
jako model predykcyjny sygnałów biologicznych do
wykrywania ich zaburzeń (np. detekcji stanu fibrylacji
komór serca implantowane defibrylatory)
69
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Uśrednianie synchroniczne sygnału
Filtry dolnoprzepustowe nie są skuteczne przy
wygładzaniu sygnałów, w których pasmo zakłóceń
pokrywa się z pasmem składowych użytecznych sygnału
x(t)= s(t)+ n(t)
x(t)= s(t)+ n(t)
A(f) A(f)
sygnał
sygnał
użyteczny
użyteczny
zakłócenie zakłócenie
S(�)
S(�)
N(�)
N(�)
f f
Widma częstotliwościowe sygnału i zakłócenia
70
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Uśrednianie synchroniczne sygnału
Idea uśredniania synchronicznego
Sygnały synchronizujące
71
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Uśrednianie synchroniczne sygnału
Uśrednianie synchroniczne sygnału jest skuteczne przy
spełnieniu następujących warunków:
składowe deterministyczne sygnału powinny
występować okresowo (niekoniecznie w regularnych
występować okresowo (niekoniecznie w regularnych
odstępach)
sygnał zakłócający powinien być sygnałem losowym,
nieskorelowanym ze składowymi deterministycznymi
sygnału
powinna istnieć możliwość detekcji cech sygnału
potrzebnych do synchronizacji kolejnych jego cykli
72
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Uśrednianie synchroniczne sygnału
N
1
x = xn
Ć
"
N
n=1
& &
^
x
x1 x2 x3 xN
73
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Uśrednianie synchroniczne sygnału
Odchylenie standardowe sygnału: �s
Odchylenie standardowe zakłócenia: �n
�
s
Stosunek sygnału do zakłócenia:
SNR =
�
�
n
n
Po N uśrednieniach:
2
10
�
s
SNRN = N
�
n
1
10
Zatem poprawa SNR po
N
N uśrednieniach wynosi:
0
10
74
0 100 200 300 400 500
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Uśrednianie synchroniczne sygnału
Zastosowania:
detekcja podszumowa sygnału tj. dla �s << �n
(zastosowania w telekomunikacji)
analiza elektrycznych potencjałów wywołanych
analiza elektrycznych potencjałów wywołanych
mózgu, tj. potencjałów generowanych w mózgu o
amplitudzie kilku mikrowoltów na skutek okresowego
pobudzenia bodz
cem: świetlnym (potencjały
wzrokowe), dz
więkowym (potencjały słuchowe) lub
dotykowym (potencjały czuciowe)
75
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Uśrednianie synchroniczne sygnału
Ćwiczenie:
Ćwiczenie: Zakłócony sygnał EKG  ecg_s.mat
poddaj wygładzeniu przez synchronizowane
sumowanie jego kolejnych cykli.
Sprawdz
jaką poprawę SNR uzyskałeś po
Sprawdz
jaką poprawę SNR uzyskałeś po
kolejnych 2, 3, 4, & uśrednieniach sygnału.
76
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Filtracja medianowa sygnałów
Mediana  element środkowy w ciągu
uporządkowanych liczb, np.:
x(n)={1, 5, -7, 101, -25, 3, 0, 11, 7}
x(n)={1, 5, -7, 101, -25, 3, 0, 11, 7}
Porządkowanie (sortowanie) ciągu:
xs(n)={-25, -7, 0, 1, 3, 5, 7, 11, 101 }
Element środkowy
77
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Filtracja medianowa sygnałów
sortowanie
czas dyskretny
czas dyskretny
M=2m+1=5 M=5
M=2m+1=5 M=5
78
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Filtracja medianowa sygnałów
Filtracja medianowa sygnału:
2 5 -1 7 0 3 -15 9 1 20 6 -9
xn
Kierunek przesuwania
-1 2 5 7
1 2
okna filtracji
med3(xn)
2 2 5 0 3 0 ? ? ? ? 6 6
wyznacz medianę
79
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Filtracja medianowa sygnałów
Własności filtru medianowego:
usuwają zakłócenia impulsowe sygnału o szerokości
mniejszej niż połowa szerokości okna filtracji
nie zniekształca położenia i kształtu zbocz impulsów
o szerokości większej niż połowa okna filtracji
zakres [min, max] sygnału po filtracji medianowej nie
jest przekroczony
z uwagi na konieczne sortowanie próbek sygnału w
oknie filtracji koszt obliczeniowy filtracji medianowej
jest wielokrotnie większy od filtracji splotowej
80
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Filtracja medianowa sygnałów
Ćwiczenie:
Ćwiczenie: Utwórz przebieg prostokątny
( help square ) i dodaj losowe zakłócenia
impulsowe do tego przebiegu.
impulsowe do tego przebiegu.
Porównaj wynik filtracji tego sygnału uzyskany po
liniowej filtracji wygładzającej z wynikiem filtracji
uzyskanym po zastosowaniu filtru medianowego
( help medfilt1 ).
81
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Filtracja medianowa -przykład
400
390
x(t)
380
y=med(x)
370
360
350
340
330
M=7
320
45
310
40
0 20 40 60 80 100 120 140 160
35
30
25
y=x-med(x)
20
15
10
5
82
0
-50
20 40 60 80 100 120 140 160
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Filtracja medianowa sygnałów
Ćwiczenie:
Ćwiczenie: Zastosuj ideę filtracji medianowej
pokazanej na slajdzie 31 do detekcji zespołów QRS w
sygnale  ecg_mit.mat .
sygnale  ecg_mit.mat .
Jaka jest optymalny zakres szerokości okna filtracji dla
tego zastosowania filtru medianowego?
83
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej
Filtracja nieliniowa sygnałów
Czy dla filtrów nieliniowych można wyznaczyć
transmitancje filtru?
84
Paweł Strumiłło, Instytut Elektroniki Politechniki A
ódzkiej