Laboratorium Cyfrowego Przetwarzania Sygnałów
Ćwiczenie 6 – Filtracja optymalna i adaptacyjna
- 1 -
Wydział Elektryczny
Zakład Automatyki
LABORATORIUM CYFROWEGO PRZETWARZANIA SYGNAŁÓW
Ćwiczenie 6
Filtracja optymalna i adaptacyjna
1. Cel ćwiczenia
• Zapoznanie się z projektowaniem filtrów optymalnych metodą minimalizacji wskaźnika jakości i
algorytmami LMS i RLS adaptacji współczynników filtrów.
• Przeprowadzenie eksperymentów symulacyjnych z filtrami adaptacyjnymi w środowisku
Simulink.
2. Podstawy teoretyczne
Zagadnienie estymacji (odtworzenia) użytecznego sygnału s(n) na podstawie innego dostępnego
sygnału x(n), który jest zaszumioną wersją s(n), tzn. x(n) = s(n) + w(n) (w(n) – niepożądane zakłócenie)
jest jednym z centralnych problemów cyfrowego przetwarzania sygnałów. Trudność polega na tym, że
widmo zakłócenia pokrywa się zwykle z widmem sygnału użytecznego i filtracja zakłócenia jest
związana z usunięciem części sygnału użytecznego. Dlatego filtrów optymalnych nie projektuje się w
dziedzinie częstotliwości, ale z wykorzystaniem statystycznych właściwości sygnału podlegającego
filtracji. W ćwiczeniu będziemy zajmować się estymacją liniową przy addytywnych zakłóceniach.
2.1. Filtracja optymalna
Struktura filtra optymalnego jest przedstawiona na rys.1. Ma on dwa sygnały wejściowe: sygnał
filtrowany x(n) oraz sygnał odniesienia d(n), i dwa sygnały wyjściowe: wynik filtracji y(n), który jest
estymatą sygnału odniesienia, tj.
ˆ
( )
( )
y n
d n
=
, oraz błąd estymacji e(n) = d(n) - y(n). Zakładamy, że
sygnały wejściowe są stacjonarne (w szerszym sensie). Ograniczymy się do filtrów nierekursywnych (o
skończonej odpowiedzi impulsowej SOI (FIR)):
1
2
0
1
2
( )
M
M
H z
h
h z
h z
h z
−
−
−
=
+
+
+
+
…
,
(6.1)
ponieważ gwarantują one stabilność i dlatego są najczęściej stosowane jako filtry adaptacyjne. Wyjście
filtra można zapisać w formie:
0
( )
(
)
( )
M
T
k
k
y n
h x n k
n
=
=
−
=
∑
h x
,
(6.2)
gdzie:
0
1
2
[ , , , ,
]
T
M
h h h
h
=
h
…
,
( ) [ ( ), (
1), (
2), , (
)]
T
n
x n x n
x n
x n M
=
−
−
−
x
…
(wskaźnik
T
oznacza
transpozycję).
Zadaniem filtra jest takie przekształcanie sygnału wejściowego x(n), aby wynik filtracji był jak
najbardziej zbliżony do sygnału odniesienia przy założonym kryterium błędu. Najczęściej stosuje się
kryterium minimalizacji błędu średniokwadratowego (mean-square error MSE)
2
[ ( )]
J
E e n
=
, (E[
⋅] oznacza wartość oczekiwaną)
(6.3)
gdzie
0
( )
( )
( )
( )
(
)
M
k
k
e n
d n
y n
d n
h x n k
=
=
−
=
−
−
∑
.
(6.4)
Laboratorium Cyfrowego Przetwarzania Sygnałów
Ćwiczenie 6 – Filtracja optymalna i adaptacyjna
- 2 -
Rys. 1. Struktura filtra optymalnego. W przypadku stacjonarnym parametry filtra H(z) są wyznaczane
jeden raz przez minimalizację kryterium błędu związanego z e(n)
Błąd dopasowania jest funkcją współczynników filtra {h
k
}:
(
)
2
2
( )
( )
( )
( )
( )
[ ( )] 2
T
T
n
T
n
dx
xx
J
E d n
n
E d n
⎡
⎤
=
−
=
−
+
⎢
⎥
⎣
⎦
h
h x
h r
h R h
,
(6.5)
gdzie
( )
[ ( ) ( )]
n
T
xx
E
n
n
=
R
x
x
jest macierzą autokorelacji sygnału wejściowego
( )
( )
,
(
)
[ (
) (
)],
,
0,1, ,
n
n
xx
xx
i j
R
i
j
E x n i x n j
i j
M
⎡
⎤ =
−
=
−
−
=
⎣
⎦
R
…
,
(6.6)
a
( )
[ ( ) ( )]
n
dx
E d n
n
=
r
x
- wektorem korelacji wzajemnej pomiędzy d(n) i x(n-i):
( )
( )
( )
[ ( ) (
)],
0,1, ,
n
n
xx
dx
i
r
i
E d n x n i
i
M
⎡
⎤ =
=
−
=
⎣
⎦
r
…
(6.7)
Optymalne wartości współczynników filtra otrzymuje się w wyniku przyrównania do zera pochodnej:
( )
( )
( )
2
2
T
n
T
n
dx
xx
J
∂
= −
+
=
∂
h
h r
h R h 0
h
,
(6.8)
co prowadzi do układu liniowego tzw. równań normalnych:
( )
( )
n
n
xx
dx
⋅ =
R
h r
(6.9)
Rozwiązanie tych równań daje optymalne współczynniki filtra:
( )
1 ( )
( )
[
]
, det
0
opt
n
n
n
xx
dx
xx
−
=
≠
h
R
r
R
(6.10)
1
( )
( )
( )
( )
0
( )
( )
( )
( )
1
( )
( )
( )
( )
(0)
(1)
( )
(0)
(1)
(0)
(
1)
(1)
( )
(
1)
(0)
( )
n
n
n
n
opt
xx
xx
xx
dx
n
n
n
n
opt
xx
xx
xx
dx
n
n
n
n
opt
xx
xx
xx
dx
M
R
R
R
M
r
h
R
R
R
M
r
h
R
M
R
M
R
r
M
h
−
⎡
⎤
⎡
⎤
⎡
⎤
⎢
⎥
⎢
⎥
⎢
⎥
−
⎢
⎥
⎢
⎥
⎢
⎥ =
⋅
⎢
⎥
⎢
⎥
⎢
⎥
⎢
⎥
⎢
⎥
⎢
⎥
−
⎢
⎥ ⎢
⎥
⎢
⎥
⎣
⎦ ⎣
⎦
⎣
⎦
Warunek odwracalności macierzy
( )
det
0
n
xx
≠
R
jest nazywany warunkiem wystarczającego pobudzenia.
Filtr (6.10) jest nazywany filtrem Wienera. Dla filtracji optymalnej zachodzi zależność:
( )
( )
2
( )
2 [ ( ) (
)] 0,
0,1, ,
k
k
J
e n
E e n
E e n x n k
k
M
h
h
⎡
⎤
∂
∂
=
= −
−
=
=
⎢
⎥
∂
∂
⎣
⎦
h
…
, (6.11)
co oznacza, że filtr usuwa korelację błędu estymacji ze wszystkimi poprzednimi (i bieżącą) próbkami
sygnału wejściowego:
( )
( )
[ ( ) (
)] 0
n
ex
r
k
E e n x n k
=
−
=
.
(6.12)
Jest to tzw. właściwość ortogonalności. Jej geometryczna interpretacja jest pokazana na rys.2.
H(z)
x(n)
błąd
e(n)
_
+
y(n)
sygnał odniesienia
d(n) = s(n)
zakłócenie
w(n)
Laboratorium Cyfrowego Przetwarzania Sygnałów
Ćwiczenie 6 – Filtracja optymalna i adaptacyjna
- 3 -
Rys. 2. Geometryczna interpretacja właściwości ortogonalności estymaty otrzymanej metodą błędu
średniokwadratowego MSE dla dwóch wymiarów (M=2). Wyjście filtra
ˆ
( )
( )
y n
d n
=
.
2.2. Estymacja współczynników filtra optymalnego z próbek
W zagadnieniach praktycznych teoretyczne kryterium (6.3) oraz korelacje (6.6) i (6.7) wyznaczane z
wartości oczekiwanych zastępuje się ich estymatami
ˆ
xx
R
i
ˆ
dx
r
obliczonymi z próbek sygnałów x(n) i d(n).
Załóżmy, że mamy do dyspozycji po N próbek sygnału wejściowego i sygnału odniesienia:
[ (0), (1), (2), , (
1)]
T
x
x
x
x N
=
−
x
…
,
[ (0), (1), (2), , (
1)]
T
d
d
d
d N
=
−
d
…
,
a próbki na wyjściu filtra określa równanie (6.2). Ponieważ zwykle liczba próbek jest większa niż długość
filtra N>M, otrzymuje się układ nadokreślony, tzn. liczba równań jest większa niż liczba niewiadomych
(poszukiwanych współczynników filtra). Istnieją dwie metody obliczeń w zależności od zakresu
sumowania we wzorze na estymatę błędu dopasowania. Każda z nich daje inne wartości estymat
ˆh
współczynników filtra optymalnego.
1. Metoda autokorelacji
Kryterium błędu:
1
2
0
( ),
( )
( )
( )
N M
n
n
J
e n
e n
d n
y n
+ −
=
=
=
−
∑
(6.13)
W metodzie tej mamy N+M równań na wartości wyjściowe y(n). Kolejne próbki y(n) są obliczane jako
suma iloczynów współczynników filtra i odpowiadających im próbek x(n), jak na poniższym rysunku. Jak
widać, ciąg próbek sygnału wejściowego jest na krańcach uzupełniany zerami. Poza pierwszymi N
równaniami w pozostałych po prawej stronie występują zera.
Metoda autokorelacji prowadzi do symetrycznej macierzy równań normalnych i gwarantuje, że filtr H(z)
jest minimalnofazowy.
x
0
x
1
x
2
x
n-M
x
n-M+1
x
n-1
x
n
x
N-1
0
0
0
0 0
h
M
h
1
h
0
.....
...... .....
...
....
h
M
h
1
h
0
.....
h
M
h
1
h
0
.....
....
d
0
d
1
d
2
d
n-M
d
n-M+1
d
n-1
d
n
d
N-1
0
0 0
..... ...
....
....
Laboratorium Cyfrowego Przetwarzania Sygnałów
Ćwiczenie 6 – Filtracja optymalna i adaptacyjna
- 4 -
2. Metoda kowariancji
Kryterium błędu:
1
2
( ),
( )
( )
( )
N
n
n M
J
e n
e n
d n
y n
−
=
=
=
−
∑
(6.14)
W metodzie tej mamy N-M równań na wartości wyjściowe y(n). Kolejne próbki y(n) są obliczane jako
suma iloczynów współczynników filtra i odpowiadających im próbek x(n), jak na poniższym rysunku.
Tworzenie równań jest ograniczone do takiego zakresu, że współczynniki filtra "nie wychodzą" poza
rekord próbek wejściowych. Zwróćmy uwagę, że nie wykorzystuje się informacji z pierwszych M próbek
sygnału odniesienia.
Metoda kowariancji nie gwarantuje minimalnofazowości H(z).
2.3. Filtracja adaptacyjna
Filtracja adaptacyjna polega na dodaniu do struktury filtra optymalnego algorytmu adaptacji
współczynników H(z) (rys. 3), które są funkcjami czasu n i są na bieżąco dopasowywane w taki sposób,
aby zapewniać optymalną, w sensie przyjętego kryterium błędu, filtrację sygnału wejściowego.
Rys. 3. Struktura filtra adaptacyjnego. Parametry h(n) filtra H
n
(z) są wyznaczane na bieżąco według
określonego algorytmu adaptacyjnego minimalizującego kryterium błędu związane z e(n)
Jeżeli wszystkie sygnały w układzie są stacjonarne, to po okresie adaptacji współczynniki filtra powinny
się ustalić i przyjąć wartości optymalne (warunek zbieżności procesu adaptacji):
lim ( )
opt
n
n
→∞
=
h
h
(6.15)
W adaptacyjnym przetwarzaniu sygnałów stosuje się metody optymalizacji wieloparametrycznej
(optymalizujemy M+1 parametrów transmitancji), które polegają poszukiwaniu minimów
deterministycznych wieloargumentowych funkcji kosztu J(
⋅). Najczęściej stosowane są metody
gradientowe, w których modyfikacja
Δh(n) jest w każdej chwili proporcjonalna do wektora ujemnego
gradientu funkcji kosztu:
( )
1
( )
(
1)
( )
(
1)
2
n
n
n
n
n
n
=
− + Δ
=
− − μ
⋅∇
h
h
h
h
W
(6.16)
0
1
( ( ))
( )
( )
( )
,
,
,
,
( )
( )
( )
( )
T
n
M
J
n
J n
J n
J n
n
h n
h n
h n
⎡
⎤
∂
∂
∂
∂
∇ =
= ⎢
⎥
∂
∂
∂
∂
⎣
⎦
h
h
…
(6.17)
H
n
(z)
x(n)
błąd
e(n)
_
+
y(n)
sygnał odniesienia
d(n)
Algorytm
adaptacji
x
0
x
1
x
n-M
x
n-M+1
x
n-1
x
n
x
N-M-1
h
M
h
1
h
0
.....
.... ......
...
....
h
M
h
1
h
0
.....
h
M
h
1
h
0
.....
...
x
M-1
x
M
x
M+1
x
N-M
x
N-2
x
N-1
d
0
d
1
d
n-M
d
n-M+1
d
n-1
d
n
d
N-M-1
.... ......
...
....
...
d
M-1
d
M
d
M+1
d
N-M
d
N-2
d
N-1
Laboratorium Cyfrowego Przetwarzania Sygnałów
Ćwiczenie 6 – Filtracja optymalna i adaptacyjna
- 5 -
W
(n)
jest dodatkową macierzą wagową poprawiającą zbieżność (np. odwrotność hesjanu
( )
2
1
[
]
n
n
−
= ∇
W
w procedurze Newtona-Raphsona). Współczynnik adaptacji
μ decyduje o szybkości przestrajania (Δh(n)
jest proporcjonalne do
μ). Ze względu na zbieżność algorytmy zachodzi warunek ograniczający:
max
2
0
< μ <
λ
,
(6.18)
gdzie
λ
max
jest największa wartością własną macierzy autokorelacji sygnału wejściowego
( )
n
xx
R
.
2.4. Filtry adaptacyjne LMS
W algorytmie adaptacyjnym LMS (Least Mean Squares) kryterium błędu
2
( )
n
J
e n
=
.
(6.19)
Zadaniem jest więc minimalizacja chwilowej (a nie oczekiwanej) wartości błędu kwadratowego, dlatego
filtrację nazywa się filtracją bez pamięci albo optymalizacją stochastyczną. Jednak asymptotycznie dla
n
→∞ algorytm ten minimalizuje średni błąd kwadratowy.
Metoda ta obejmuje szeroką rodzinę algorytmów opisanych ogólną zależnością:
( )
(
1)
( ) ( ) ( ) ( )
n
n
n
n e n
n
=
− + μ
h
h
W
x
(6.20)
Przypadki szczególne:
1. Filtr LMS
( )
(
1)
( ) ( )
n
n
e n
n
=
− + μ
h
h
x
(6.21)
( )
( )
( )
e n
d n
y n
=
−
( )
(
1) ( )
T
y n
n
n
=
−
h
x
Zaletą algorytmu jest prostota i mała złożoność obliczeniowa, a wadą – wolna zbieżność algorytmu.
2. Unormowany filtr LMS (NLMS)
( )
(
1)
( ) ( ) ( )
n
n
n e n
n
=
− + μ
h
h
x
,
(6.22)
gdzie
2
0
( )
( ) ( )
(
)
M
T
k
n
n
n
x n k
=
μ
μ
μ
=
=
γ +
γ +
−
∑
x
x
(mały parametr
γ zapobiega zerowaniu się mianownika).
Filtry NLMS charakteryzują się szybszą zbieżnością i lepszą stabilnością w porównaniu z LMS.
3. Zdekorelowany filtr LMS
Zamiast x(n) na wejście podawany jest zdekorelowany sygnał
( ) (
1)
( )
( )
( ) (
1),
( )
(
1) (
1)
T
T
n
n
n
n
n
n
n
n
n
−
=
− α
−
α
=
−
−
x
x
v
x
x
x
x
(6.23)
który ma stosunek
λ
min
/
λ
max
bliski jedności, co zwiększa szybkość adaptacji.
α(n) jest współczynnikiem
korelacji wektorów x(n) i x(n-1).
2.5. Filtry adaptacyjne RLS
Popularnym algorytmem adaptacji jest algorytm LS (Least Squares) z kryterium błędu w postaci sumy
kwadratów:
2
0
( )
n
n
k
J
e n
=
=
∑
(6.24)
lub ważonej sumy kwadratów WLS (Weighted LS):
2
0
( )
n
n k
n
k
J
e n
−
=
=
λ
∑
,
(6.25)
gdzie
λ jest współczynnikiem zapominania starych błędów o wartości z zakresu
0.95
1
≤ λ ≤
. Algorytm
z wagami umożliwia adaptację do zmieniających się charakterystyk sygnałów, kiedy dane z przeszłości
przestają być aktualne. Ponieważ kryterium zależy od sumy błędów, filtry LS nazywa się filtrami z
Laboratorium Cyfrowego Przetwarzania Sygnałów
Ćwiczenie 6 – Filtracja optymalna i adaptacyjna
- 6 -
pamięcią. W trakcie adaptacji ruch trajektorii układu w przestrzeni parametrów w kierunku minimum jest
mniej chaotyczny niż w przypadku metody LMS.
Równanie zmian współczynników można zapisać w postaci:
( )
(
1)
( ) ( )
n
n
n e n
=
− +
h
h
K
,
(6.26)
gdzie
( )
1
ˆ
( ) [
]
( )
n
xx
n
n
−
=
K
R
x
nazywa się wzmocnieniem Kalmana.
Ze względu na konieczność odwracania w każdym kroku macierzy autokorelacji algorytm LS
implementuje się w formie rekurencyjnej RLS (Recursive LS)
( )
(
1)
( ) ( )
n
n
n e n
=
− +
h
h
K
( )
( )
( )
e n
d n
y n
=
−
,
( )
(
1) ( )
T
y n
n
n
=
−
h
x
1
1
(
1) ( )
( )
( ) ( )
1
( ) (
1) ( )
T
n
n
n
n
n
n
n
n
−
−
λ
−
=
=
+ λ
−
P
x
K
P
x
x
P
x
, (gdzie
( )
1
ˆ
( ) [
]
n
xx
n
−
=
P
R
)
(6.27)
1
1
( )
[
( ) ( )] (
1)
T
n
n
n
n
−
−
= λ
− λ
−
P
I
K
x
P
, (6.28)
z warunkami początkowymi:
0
2
1
(0)
,
(0)
ˆ
x
=
=
σ
P
I
h
h
.
Algorytm RLS jest bardziej złożony obliczeniowo niż LMS, ale zwykle daje szybszą i płynniejszą
zbieżność.
2.6. Typowe zastosowania filtracji adaptacyjnej
Poniżej przedstawiono kilka typowych struktur zastosowania filtracji adaptacyjnej ilustrujących
interpretację poszczególnych sygnałów.
We wszystkich przypadkach wyjście filtra y(n) jest estymatą tylko tej składowej sygnału odniesienia
d(n), która jest skorelowana z sygnałem wejściowym x(n). Jeżeli d(n) = s(n) + x
1
(n), gdzie tylko x
1
(n) jest
skorelowane z x(n), to wyjście
1
ˆ
( )
( )
y n
x n
=
.
Rys. 4. Identyfikacja adaptacyjna. Transmitancja filtra adaptacyjnego jest estymatą transmitancji
identyfikowanego układu
ˆ
( )
( )
H z
G z
=
.
Rys. 5. Adaptacyjna predykcja sygnału
H
n
(z)
x(n)=d(n-k)
błąd predykcji
e(n)
_
+
ˆ
( )
( )
y n
d n
=
sygnał odniesienia
d(n) = s(n)
Predyktor k-krokowy
z
-k
ˆ
( )
( )
y n
d n
=
H
n
(z)
x(n)
e(n)
_
+
Filtr adaptacyjny
- model układu
G(z)
Nieznany układ
zakłócenie w(n)
d(n)
Laboratorium Cyfrowego Przetwarzania Sygnałów
Ćwiczenie 6 – Filtracja optymalna i adaptacyjna
- 7 -
Rys. 6. Adaptacyjne kasowanie szumu (usuwanie korelacji)
Rys. 6. Przykład zastosowania adaptacyjnego kasowania szumu. Zakłócenie skorelowane z szumem
własnym okrętu dociera do sonaru jako jego stłumiona i opóźniona wersja w(n)=
α⋅w
0
(n-n
0
). Zadaniem
filtracji adaptacyjnej jest wytłumienie tej składowej w sygnale sonaru.
Literatura
1. Stranneby D.: „Cyfrowe przetwarzanie sygnałów. Metody, algorytmy, zastosowania”, Wyd. BTC, 2004.
2. Zieliński T.P.: „Cyfrowe przetwarzanie sygnałów. Od teorii do zastosowań”, WKŁ, 2005.
'
#
szum własny w
0
(n)
sygnał sonaru d(n)=s(n)+w(n)
echo obiektu s(n)
H
n
(z)
x(n)=w
0
(n)
ˆ
( )
( )
e n
s n
=
_
+
d(n)=s(n)+w(n)
Filtr adaptacyjny
Źródło
szumu
Źródło
sygnału
ˆ
( )
( )
y n
w n
=
Laboratorium Cyfrowego Przetwarzania Sygnałów
Ćwiczenie 6 – Filtracja optymalna i adaptacyjna
- 8 -
3. Obliczenia komputerowe - zadania do wykonania
Do obliczeń wykorzystywana jest biblioteka symulacyjna Simulinka DSP Blockset oraz biblioteka
Matlaba Signal Processing Toolbox.
3.1. Liniowa predykcja adaptacyjna
A. Otworzyć model symulacyjny .mdl do predykcji liniowej:
>> c6_predict
Liniowa predykcja adaptacyjna
Osc
Delay
Flip
ZOH
User
Wsp. predykcji h(n)
t
Wektor chwil czasu
Szum pomiarowy
Signal
Generator
x
S.wejsciowy
d
S. odniesienia
dp
Predykcja
In
Err
Out
Taps
nLMS
LMS
Adaptive Filter
z
-2
Clock
d(n)
y(n) = dp(n)
e(n)
x(n)
• Ustawić (sprawdzić) parametry symulacji w menu Simulation | Parameters:
Stop time=1000, Solver options: Fixed-step, discrete
Ustawić (sprawdzić) parametry bloków:
Signal generator: Waveform: sine, Amplitude=1, Frequency=0.015 Hz,
Sample time=1 (okres próbkowania)
Szum: Variance=0.0 (brak zakłócenia), Sample time=1
Delay=2 (opóźnienie sygnału wejściowego filtra w próbkach, predykcja 2-krokowa)
LMS Filter: FIR length=16, Step-size mu=1.0, Use normalization: Tak
(długość filtra L, współczynnik adaptacji
μ, algorytm NLMS z normalizacją)
• Otworzyć okno oscyloskopu OSC i przeprowadzić symulację (przycisk Play lub Simulation |
Start). Obserwować przebieg adaptacji.
• Zapisać wykres z końcową charakterystyką impulsową filtra predykcyjnego LMS h(n) (Copy Figure).
W oknie komend uruchomić program skryptowy do wykreślania sygnałów zarejestrowanych w
pamięci:
» c6_predictplot(t,d,x,dp)
Wkreślane są: sygnał odniesienia d(n), wejście filtra x(n) (opoźniona i zaszumiona wersja sygnału
odniesienia), wyjście filtra
ˆ
( )
( )
y n
d n
=
(predykcja) oraz błąd predykcji
ˆ
( )
( )
( )
e n
d n
d n
=
−
. Zapisać okno
z wykresami.
B. Powtórzyć obliczenia i zarejestrować wyniki dla pozostałych 3 kombinacji parametrów: wariancja
szumu zakłócenia Variance=0.0, 0.05, opóźnienie Delay=2, 10.
C. Porównać przebiegi adaptacji dla jednej z kombinacji parametrów z szumem dla innej wartości
współczynnika adaptacji, np. mu=0.3 (filtr LMS). Sprawdzić, jak działa adaptacja LMS bez
normalizacji (usunąć zaznaczenie opcji Use normalization w filtrze).
&
Przeanalizować uzyskane wyniki i wyciągnąć wnioski.
3.2. Estymacja opóźnienia
A. Otworzyć model symulacyjny .mdl do predykcji liniowej:
Laboratorium Cyfrowego Przetwarzania Sygnałów
Ćwiczenie 6 – Filtracja optymalna i adaptacyjna
- 9 -
>> c6_delay
Estymacja opoznienia
Flip
y
Wyjscie y(n)
User
Wsp.filtra
t
Wektor chwil czasu
x
Wejscie x(n)
Szum pomiarowy
Sine Wave
d
S.odniesienia d(n)
Results
In
Err
Out
Taps
nLMS
LMS
Adaptive Filter
1
Kx
0.5
z
-4
Delay
Clock
Echo d(n)
Error e(n)
Input x(n)
y(n)
• Ustawić (sprawdzić) parametry symulacji w menu Simulation | Parameters:
Stop time=0.5, Solver options: Fixed-step, discrete
Ustawić (sprawdzić) parametry bloków:
Sine wave: Amplitude=1, Frequency=50 Hz, Sample time=1e-3 (okres
próbkowania T
s
=1 ms)
Szum: Variance=0.2, Sample time=1e-3
Delay=4 (opóźnienie echa – składowej sygnału odniesienia)
LMS Filter: FIR length=16, Step-size mu=1.0, Use normalization: Tak
(długość filtra L, współczynnik adaptacji
μ, algorytm NLMS z normalizacją)
• Otworzyć okno oscyloskopu OSC i przeprowadzić symulację (przycisk Play lub Simulation |
Start). Obserwować przebieg adaptacji.
• Zapisać wykres z końcową charakterystyką impulsową filtra predykcyjnego LMS h(n) - estymatę
opóźnienia (Copy Figure). W oknie komend uruchomić program skryptowy do wykreślania sygnałów
zarejestrowanych w pamięci:
» c6_delayplot(t,d,x,y)
Wkreślane są: sygnał odniesienia d(n) (wejście filtra + echo), wejście filtra x(n) (zaszumiony sygnał
harmoniczny), wyjście filtra y(n), błąd estymacji
( )
( )
( )
e n
d n
y n
=
−
. Zapisać okno z wykresami.
B. Powtórzyć obliczenia i zarejestrować wyniki dla czasu opóźnienia Delay=10.
C. Porównać przebieg adaptacji dla innej wartości współczynnika adaptacji, np. mu=0.3 (filtr LMS).
D. Powtórzyć obliczenia i zarejestrować wyniki dla samego echa jako sygnału odniesienia : parametr
Kx=0 (przerwanie gałęzi przepływu sygnału).
&
Zwrócić uwagę na estymatę opoźnienia h(n) (parametry filtra). Wyjaśnić wynik.
E. Przywrócić Kx=1. Powtórzyć obliczenia i zarejestrować wyniki przy braku szumu zakłócenia
Szum: Variance=0.
&
Zwrócić uwagę na estymatę opoźnienia h(n) (parametry filtra). Wyjaśnić wynik.
3.3. Identyfikacja układu
A. Otworzyć model symulacyjny .mdl do predykcji liniowej:
>> c6_ident
Laboratorium Cyfrowego Przetwarzania Sygnałów
Ćwiczenie 6 – Filtracja optymalna i adaptacyjna
- 10 -
Identyfikacja adaptacyjna
x
wymuszenie + szum
ZOH
User
Wsp. estymatora he(n)
t
Wektor chwil czasu
fir1
Uklad identyfikowany (FIR)
Szum pomiarowy
Signal
Generator
d
S.wyjsciowy
In
Err
Out
Taps
nLMS
LMS
Adaptive Filter
FFT
Freq Response
y
Estymacja y
e
Error e(n)
Clock
y(n) = de(n)
d(n)
x(n)
e(n)
• Ustawić (sprawdzić) parametry symulacji w menu Simulation | Parameters:
Stop time=1500, Solver options: Fixed-step, discrete
Ustawić (sprawdzić) parametry bloków:
Signal generator: square (wymuszenie prostokątne), Amplitude=1,
Frequency=0.02 Hz, Sample time=1
Szum: Variance=0.0 (brak zakłócenia), Sample time=1
Układ identyfikowany (FIR): Bandstop, Filter order=32 (filtr
pasmowozaporowy), Lower cutoff freq=0.3, Upper cutoff freq=0.7, Hamming
window (częstotliwość Nyquista = 1).
LMS Filter: FIR length=32, Step-size mu=1.0, Use normalization: Tak
(długość filtra L, współczynnik adaptacji
μ, algorytm NLMS z normalizacją)
• Otworzyć okno oscyloskopu OSC i przeprowadzić symulację (przycisk Play lub Simulation |
Start). Obserwować przebieg adaptacji.
• Zapisać wykres z końcową charakterystyką impulsową filtra LMS h(n) - model dynamiki
identyfikowanego układu oraz jego charakterystykę amplitudową (Copy Figure). (Uwaga: w oknie
charakterystyki filtra LMS częstotliwość Nyquista = 0.5 – połowa częstotliwości próbkowania). W
oknie komend uruchomić program skryptowy do wykreślania sygnałów zarejestrowanych w pamięci:
» c6_identplot(t,d,x,y)
B. Powtórzyć obliczenia i zarejestrować wyniki z szumem zakłócenia: Szum: Variance=0.2.
C. Zmienić identyfikowany układ na pasmowoprzepustowy – wybrać z listy Filter type:
Bandpass. Przeprowadzić symulację (z szumem) i zarejestrować wyniki.
&
Czy charakterystyka częstotliwościowa filtra adaptacyjnego we właściwy sposób odwzorowuje
charakterystykę identyfikowanego układu?
D. Ustawić długi czas symulacji Simulation | Parameters: Stop time=50000.
Uruchomić symulację i obserwować przebieg adaptacji przy zmienianiu parametrów on-line. Można
zmieniać typ (np. BP
→BS→LP) i częstotliwości graniczne układu, wariancję szumu, współczynnik
adaptacji mu (np. mu=0.3) filtra adaptacyjnego, sygnał wymuszenia (wypróbować np. wymuszenie
sinusoidalne sine, zamiast prostokątnego square, przy braku szumu, tj. Variance=0.0).
&
Przeanalizować uzyskane wyniki i wyciągnąć wnioski.
E. Zamienić filtr adaptacyjny w schemacie blokowym na filtr RLS: zaznaczyć blok filtra LMS i usunąć
go, otworzyć bibliotekę bloków Simulinka Library Browser (przycisk na belce okna schematu) i wybrać
Laboratorium Cyfrowego Przetwarzania Sygnałów
Ćwiczenie 6 – Filtracja optymalna i adaptacyjna
- 11 -
DSP Blockset | Filtering | Adaptive filters | RLS Adaptive filter i przeciągnąć go w miejsce usuniętego
bloku (Uwaga: strzałki wejść i wyjść muszą się połączyć z portami bloku!). Parametry filtra RLS:
RLS Filter: FIR length=32, Memory weighting factor=1.0 (współczynnik
zapominania
λ), Initial value of taps=0.0, Initial variance estimate=0.1
Przywrócić wartości pozostałych parametrów jak w punkcie A. Powtórzyć obliczenia jak w punktach A,
B, C i zarejestrować wyniki.
&
Porównać wyniki z uzyskanymi dla filtra NLMS.
F. (*) Powtórzyć symulację z filtrem RLS z długim czasem symulacji jak w punkcie D. Zaobserwować
efekty po zmianie Memory weighting factor filtra RLS np. na 0.99.
&
Przeanalizować uzyskane wyniki i wyciągnąć wnioski.
3.4. Adaptacyjne kasowanie szumu
A. Otworzyć model symulacyjny .mdl do kasowania szumu:
>> c6_noisecanc
Adaptacyjne kasowanie szumu - filtr LMS/RLS
y(n)
h(n)
ZOH
User
Wsp. filtra
t
Wektor chwil czasu
Szum pierwotny
Switch
Signal Gen.
d
S.odniesienia (zaszumiony)
s
S. oryginalny
y
Output y(n)
Osc
In
Err
Out
Taps
nLMS
LMS
Adaptive Filter
[ts,s]
From
Workspace
FFT
Freq Response
Flip
fir1
Filtr LP szumu
e
Estymata se(n)
Clock
s(n)
d(n) = s(n) + w(n)
e(n) = se(n)
x(n) = w(n)
• Ustawić (sprawdzić) parametry symulacji w menu Simulation | Parameters:
Stop time=500, Solver options: Fixed-step, discrete
Ustawić (sprawdzić) parametry bloków:
Signal generator: sine (wymuszenie harmoniczne),
Amplitude=1,
Frequency=0.2 Hz, Sample time=1
Switch: przełącznik w pozycji Signal Gen. (dwukrotne kliknięcie zmienia położenie styku)
Szum: Variance=1.0, Sample time=1 (Uwaga: wariancja równa 0 oznacza zerowy
sygnał wejściowy filtra adaptacyjnego)
Filtr LP szumu: Lowpass, Filter order=32, Cutoff freq=0.5, Hamming
window (częstotliwość Nyquista = 1).
LMS Filter: FIR length=32, Step-size mu=0.3, Use normalization: Tak
(długość filtra L, współczynnik adaptacji
μ, algorytm NLMS z normalizacją)
W oknie komend zainicjalizować zmienne, np.:
>> s=0; ts=0;
• Otworzyć okno oscyloskopu OSC i przeprowadzić symulację (przycisk Play lub Simulation |
Start). Obserwować przebieg adaptacji.
• Zapisać wykres z końcową charakterystyką impulsową filtra LMS h(n) - model dynamiki
identyfikowanego układu oraz jego charakterystykę amplitudową (Copy Figure). (Uwaga: w oknie
Laboratorium Cyfrowego Przetwarzania Sygnałów
Ćwiczenie 6 – Filtracja optymalna i adaptacyjna
- 12 -
charakterystyki filtra LMS częstotliwość Nyquista = 0.5 – połowa częstotliwości próbkowania). W
oknie komend uruchomić program skryptowy do wykreślania sygnałów zarejestrowanych w pamięci:
» c6_noisecancplot(t,d,s,y)
B. Rozszerzyć pasmo przepustowe filtra LP szumu zadając Cutoff freq=0.8 i powtórzyć
symulację.
&
Jak charakterystyka amplitudowa filtra adaptacyjnego ma się do widma szumu w(n) (pamiętając, że
widmo mocy P
ww
(
Ω) jest proporcjonalne do |H
LP
(e
j
Ω
)|
2
filtra szumu i że wyjście filtra ma
kompensować składową od zakłócenia w(n) w sygnale odniesienia d(n))?
C. Powtórzyć obliczenia i zarejestrować wyniki z szumem wejściowym o większej wariancji
Variance=10.0.
D. Przeprowadzić obliczenia i zarejestrować wyniki dla współczynnika adaptacji filtra LMS mu=1.0.
E. Przeprowadzić symulację i zaobserwować efekty dla filtracji LMS bez normalizacji (skasować
zaznaczenie Use normalization w bloku filtra LMS).
&
Porównać otrzymane wyniki i wyciągnąć wnioski. Jaki jest wpływ normalizacji w algorytmie LMS?
F. Przeprowadzić eksperyment z kasowaniem szumu w sygnale dźwiękowym.
• Wczytać sygnał dźwiękowy z pliku .wav o podanej nazwie, np. (*** potrzebne głosniki ***)
>> [s,ts,fs]=c6_sound('fire8'); % sygnał, czas, częstotliwość próbkowania
(sygnał jest wysyłany na urządzenie wyjściowe i rysowany jest przebieg czasowy).
W oknie komend podawana jest liczba próbek sygnału, dla fire8.wav jest to n=6588. W zależności
od długości sygnału ustawić czas symulacji w Simulation | Parameters, np. Stop time=6500.
• Przełączyć Switch w położenie From Workspace (dwukrotne kliknięcie w obszarze bloku).
Ustawić wariancję szumu Variance=3.0. Przeprowadzić symulację z wektorem sygnału x.
• Odsłuchać wyniki (*** potrzebne głosniki ***):
>> sound(d,fs); % sygnał zaszumiony
>> sound(e,fs); % błąd estymacji, estymata
ˆ( )
s n
sygnału oryginalnego po usunięciu szumu
• Zrobic i zarejestrować wykres estymaty
ˆ( )
s n
oraz różnicy
ˆ
( )
( )
s n
s n
−
:
» c6_estymplot(t,s,e);
G. Zamienić filtr adaptacyjny w schemacie blokowym na filtr RLS jak w punkcie 3.3.E (takie same
parametry filtra RLS). Powtórzyć obliczenia jak w punkcie D.
• Zmienić współczynnik zapominania λ w bloku filtra RLS Memory weighting factor z 1.0 na
np. 0.99. Przeprowadzić obliczenia i zaobserwować utratę stabilności algorytmu adaptacji.
H. (*) Przywrócić
λ=1, przełączyć Switch w położenie Signal gen. skrócić czas symulacji
(Stop time=500) i przeprowadzić symulacje z filtrem RLS z parametrami jak w punktach A,B,C.
&
Porównać wyniki uzyskane dla filtracji NLMS i RLS.
4. Opracowanie sprawozdania
W sprawozdaniu należy zawrzeć zarejestrowane wyniki eksperymentów numerycznych z odpowiednimi
opisami oraz wyjaśnieniami wskazanych problemów.