Filtry

AS/ Funkcja przejścia i filtry

Cyfrowe filtry liniowe niezmiennicze w czasie [math]x[n] \longrightarrow \boxed{LTI} \longrightarrow y[n][/math] opisuje liniowe równanie różnicowe o stałych współczynnikach

[math] \displaystyle y[n] = \sum_{k=1}^K a_k y[n-k] + \sum_{l=0}^L b_l x[n-l] [/math]

gdzie [math]a_i[/math] i [math]b_i[/math] to współczynniki, [math]x[n][/math] to sygnał wejściowy, a [math]y[n][/math] — wyjście; w ogólniejszej postaci można je zapisać jako

[math] \displaystyle \sum_{k=0}^K a_k y[n-k] = \sum_{l=0}^L b_l x[n-l] [/math]

Funkcja przejścia (transfer function)

Zastosujmy do obu stron powyższego równania transformatę [math]\mathcal{Z}[/math]:

[math] \displaystyle \mathcal{Z}\left\{\sum_{k=0}^K a_k y[n-k] \right\} = \mathcal{Z}\left\{ \sum_{l=0}^L b_l x[n-l] \right\} [/math]

[math] \displaystyle \sum_{k=0}^K a_k \mathcal{Z}\left\{ y[n-k]\right\} = \sum_{l=0}^L b_l \mathcal{Z} \left\{x[n-l]\right\} [/math]

[math] \displaystyle \sum_{k=0}^K a_k z^{-k} Y(z) = \sum_{l=0}^L b_l z^{-l} X(z) [/math]

[math] \displaystyle Y(z) \sum_{k=0}^K a_k z^{-k} = X(z) \sum_{l=0}^L b_l z^{-l} [/math]

Dla systemu przyczynowego dostajemy:

[math] \displaystyle \frac{Y(z)}{X(z)} \equiv H(z) = \frac{\sum_{l=0}^L b_l z^{-l}}{\sum_{k=0}^K a_k z^{-k}} [/math]

[math]H(z)[/math] — funkcja systemu (system function) pozwala spójnie przedstawić działanie filtra LTI na sygnał [math]x[/math] w przestrzeni transformaty [math]\mathcal{Z}[/math]:

[math] \displaystyle Y[z]=H[z]X[z]=\frac{b_0 + b_1 z^{-1}+\dots +b_{L} z^{-L}}{a_0+a_1 z^{-1}+\dots +a_{K} z^{-K}}X[z] [/math]

W dziedzienie czasu mamy

[math] \displaystyle \begin{array}{rl} y[n] = &b_0 x[n]+ b_1 x[n-1] + \dots + b_L x[n-L]\\ &- a_1 y[n-1] - \dots - a_K y[n-K] \end{array} [/math]

Finite Impulse Response (FIR) — filtr o skończonej odpowiedzi impulsowej

Jeśli w równaniu

[math] \displaystyle \sum_{k=0}^K a_k y[n-k] = \sum_{l=0}^L b_l x[n-l] [/math]

położymy [math]a_i = 0[/math] poza [math]a_0=1[/math], dostaniemy

[math] \displaystyle y[n] = \sum_{l=0}^L b_l x[n-l] [/math]

W funkcji przejścia mianownik będzie stały i dostaniemy

[math] \displaystyle Y[z]=H[z]X[z]=\left(b_0+b_1 z^{-1}+\dots +b_L z^{-L}\right) X[z] [/math]

a w dziedzienie czasu

[math] \displaystyle y(n) = b_0 *x[n] + b_1 *x[n-1] + \dots + b_L *x[n-L] [/math]

jeśli współczynniki [math]b_i[/math] zapiszemy jako [math]b[i][/math], dostaniemy splot

[math] \displaystyle y(n) = b[0]*x[n] + b[1]*x[n-1] + \dots + b[L]*x[n-L] = b[n]*x[n] [/math]

Taki filtr nazywamy filtrem o skończonej odpowiedzi impulsowej (Finite Impulse Response, FIR), bo odpowiedź na impulsowe wzbudzenie kończy się po [math]n_b[/math] próbkach. Inna nazwa to średnia biegnąca (Moving Average, MA).

Infinite Impulse Response (IIR) — filtr o nieskończonej odpowiedzi impulsowej

Jeśli w równaniu

[math] \displaystyle \sum_{k=0}^K a_k y[n-k] = \sum_{l=0}^L b_l x[n-l] [/math]

położymy [math]b_i = 0[/math] poza [math]b_0=1[/math], dostaniemy

[math] \displaystyle \sum_{k=0}^K a_k y[n-k] = x[n] [/math]

w funkcji przejścia licznik będzie stały

[math] \displaystyle Y[z]=H[z]X[z]=\frac{1}{a_0+a_1 z^{-1}+\dots +a_{K} z^{-K}}X[z] [/math]

Operacja splotu działa tu na sekwencji wyjściowej:

[math] y[n] = x[n] - a[1]*y[n-1] - \dots - a[n_a]*y[n-n_a] [/math]

Taki filtr ten nazwać można filtrem rekursywnym lub autoregresyjnym (AR).

W praktyce filtry IIR są zwykle implementowane jako połączenie członów AR i MA, czyli:

[math] \begin{array}{ll} y[n] = b[0]*x[n] &+ b[1]*x[n-1] + \dots + b[n_b]*x[n-n_b]\\ &- a[1]*y[n-1] - \dots - a[n_a]*y[n-n_a] \end{array} [/math]

Tę wersję filtru nazywamy filtrem o nieskończonej odpowiedzi impulsowej (Infinite Impulse Response IIR) bo potencjalnie raz wzbudzony może dowolnie długo produkować niezerowe wyjście. Faza filtrowanego sygnału zaburzana jest nieliniowo (nonlinear phase filter)

Liniowe przesunięcie fazy

Jeśli współczynniki tworzą sekwencję symetryczną bądź antysymetryczną, oparty na nich filtr FIR będzie liniowo przesuwał fazę filtrowanego sygnału (linear phase filter) — sygnał filtrowany jest przesunięty w czasie o ok. [math]n_b / 2[/math].

kod (uproszczony)

import matplotlib.pyplot as plt
import numpy as np

t = np.arange(0.0, 3, 0.01)
s1 = np.sin(2 * np.pi * t)
s2= np.sin(2 * np.pi * t + np.pi/2)
s3= np.sin(2*2 * np.pi * t)
s4= np.sin(2*2 * np.pi * t + np.pi/2)
s5= np.sin(2*2 * np.pi * t + 2*np.pi/2)

fig1, ax1 = plt.subplots(3, 1)
ax1[0].plot(t, s1)
ax1[1].plot(t, s3)
ax1[2].plot(t, s1+s3, 'r')
ax1[0].set_title('sin(wt)')
ax1[1].set_title('sin(2wt)')
ax1[2].set_title('sin(wt) + sin(2wt)')

fig2, ax2 = plt.subplots(3, 1)
ax2[0].plot(t, s2)
ax2[1].plot(t, s4)
ax2[2].plot(t, s2+s4, 'r')
ax2[0].set_title('sin(wt+pi/2)')
ax2[1].set_title('sin(2wt+pi/2)')
ax2[2].set_title('sin(wt+pi/2) + sin(2wt+pi/2)')

fig3, ax3 = plt.subplots(3, 1)
ax3[0].plot(t, s2)
ax3[1].plot(t, s5)
ax3[2].plot(t, s2+s5, 'r')
ax3[0].set_title('sin(w(t+pi/2))')
ax3[1].set_title('sin(2wt+2*pi/2)')
ax3[2].set_title('sin(wt+pi/2) + sin(2wt+2*pi/2)')

plt.show()

Rzędem filtru nazywamy maksymane opóźnienie w próbkach potrzebne do wytworzenia nowej próbki wyjściowej. Dla filtrów FIR jest on równy liczbie [math]n_b[/math]. Dla filtrów IIR jest to większa z liczb [math]n_a, n_b[/math] (w praktyce często [math]n_a=n_b[/math]).

Stosując transformatę [math]Z[/math] możemy równanie z dziedziny czasu przenieść do dziedziny częstości. Filtrowanie odpowiada przemnożeniu transformaty sygnału przez transformatę funkcji przenoszenia filtru:

[math] \displaystyle Y[z]=H[z]X[z]=\frac{b[0]+b[1]z^{-1}+\dots +b[n_b]z^{-n_b}}{a[0]+a[1]z^{-1}+\dots +a[n_a]z^{-n_a}}X[z][/math]

Występująca tu funkcja [math]H(z)[/math] nosi nazwę transmitancji lub funkcja przenoszenia. Znając funkcję [math]H[/math] łatwo możemy przewidzieć co się stanie z widmem sygnału po przefiltrowaniu. Weźmy [math] z = e^{i\omega}[/math]. Wówczas transmitancja jest funkcją częstości [math]\omega[/math]. Konkretnej częstości [math]\omega_k[/math] przypisuje ona liczbę zespoloną, którą można wyrazić jako [math]A_k e^{i \phi_k}[/math].

[math] \displaystyle Y[e^{i\omega}]=H[e^{i\omega}]X[e^{i\omega}]=\frac{b[0]+b[1]e^{-i\omega}+\dots +b[n_b]e^{-i\omega n_b}}{a[0]+a[1]e^{-i\omega}+\dots +a[n_a]e^{-i\omega n_a}}X[e^{i\omega}] [/math]

W dziedzinie częstości sygnał wyrażony jest przez współczynniki Fourierowskie. Dla konkretnej częstości współczynnik taki [math]X_k = |X_k| e^{i \theta_k}[/math] (liczba zespolona) mówi z jaką amplitudą i jaką fazą exponens zespolony o danej częstości ([math]z_k = e^{i\omega_k}[/math]) wchodzi w skład sygnału.

Działanie filtru na sygnał w dziedzinie częstości polega na przemnożeniu składowej sygnału o częstości [math]f_k[/math] przez liczbę [math]A_k e^{i \phi_k}[/math]:

[math]Y(f_k) = A_k e^{i \phi_k} |X_k| e^{i \theta_k} = A_k |X_k| e^{i ( \phi_k +\theta_k)} e^{i 2\pi f_k} [/math]

W wyniku filtrowania sinusoidalna składowa sygnału o danej częstości może zmienić amplitudę i fazę ale co warto zauważyć nie zmienia częstości. Zera i bieguny filtra to odpowiednio miejsca zerowe licznika i mianownika funkcji przenoszenia.