Teoria ed elaborazione dei segnali/DFT e FFT

Da WikiAppunti.
Nota disambigua.svg Il titolo di questa pagina non è corretto per via delle caratteristiche del software MediaWiki. Il titolo corretto è B4. DFT e FFT.
Gtk-go-back-ltr.svg  B3. Trasformata di Fourier a tempo discreto Gtk-home.svg  Teoria ed elaborazione dei segnali (Elaborazione numerica dei segnali) B5. Analisi in frequenza di segnali a tempo continuo campionati  Gtk-go-forward-ltr.svg

Trasformata di Fourier discreta (DFT)[modifica]

Definizione di DFT[modifica]

4 Data una sequenza x \left( n \right) costituita da un numero N finito di campioni, la sua trasformata di Fourier discreta (DFT) X \left( k \right), all'interno del suo periodo N, è definita:

X \left( k \right) = \sum_{n = 0}^{N- 1} x \left( n \right) e^{-j2 \pi n \frac{k}{N}} , \quad  k = 0,1,2, \ldots , N-1
DFT interpretation.png

9 e può essere interpretata come la discretizzazione in frequenza della DTFT X \left( e^{j2 \pi f} \right), di cui vengono prese N frequenze f_k equispaziate:

f_k = {k \over N} , \quad  k = 0, 1, 2, \ldots , N - 1

3-7 Dal punto di vista algoritmico, la DFT è di complessità inferiore rispetto alla DTFT, e inoltre è definita in modo discreto → è rappresentabile su un calcolatore.

Inversione della DFT (IDFT)[modifica]

5 L'antitrasformata della DFT, all'interno del suo periodo N, è definita:

x \left( n \right) = \frac{1}{N} \sum_{k = 0}^{N - 1} X \left( k \right) e^{j2 \pi n \frac{k}{N}} , \quad  n = 0,1,2, \ldots, N - 1

Estensioni periodiche[modifica]

5-8 La DFT e la IDFT[1] sono periodiche di periodo N:

\begin{cases} X \left( k \right) = X \left( k + N \right) \quad k = 0 , \ldots, N-1 \\
x \left( n \right) = x \left( n + N \right) \quad n = 0 , \ldots , N-1 \end{cases}

Si definiscono le estensioni periodiche \overline{X} \left( k \right) e \overline{x} \left( n \right), rispettivamente di X \left( k \right) e x \left( n \right):

\begin{cases} \overline{X} \left( k \right) = \sum_{n = 0}^{N - 1} x \left( n \right) e^{-j2 \pi n \frac{k}{N}} \quad \forall k \\
\overline{x} \left( n \right) = \frac{1}{N} \sum_{k = 0}^{N-1} X \left( k \right) e^{j2 \pi n \frac{k}{N}} \quad \forall n \end{cases}

Tecnica di aggiunta degli zeri[modifica]

11 È possibile teoricamente ricavare per interpolazione la DTFT partendo dagli N campioni della DFT:

X \left( e^{j 2 \pi f} \right)  = \frac{1}{N} \sum_{k = 0 }^{N-1} X \left( k \right)  \frac{\sin{\left( \pi N f - k \pi \right)}}{\sin{\left( \pi f - \pi \frac{k}{N} \right)}} e^{-j \pi \left( N - 1 \right) \left( f - \frac{k}{N} \right)}
10 Dimostrazione
X \left( e^{j2 \pi f} \right) = \sum_{n = - \infty }^{+ \infty} x \left( n \right)  e^{-j2 \pi f n} =

Siccome x \left( n \right) ha per ipotesi durata finita N:

= \sum_{n = 0}^{N- 1} x \left( n \right)  e^{-j2 \pi f n} = \sum_{n = 0}^{N - 1} \left( \frac{1}{N} \sum_{k = 0}^{N-1} X \left( k  \right) \sum_{k = 0}^{N-1} X \left( k \right) e^{j2 \pi n \frac{k}{N}} \right) e^{-j2 \pi fn} = \frac{1}{N} \sum_{k = 0}^{N -1} X \left( k \right) \left( \sum_{n = 0}^{N - 1} e^{-j2 \pi f n} e^{j2 \pi n \frac{k}{N}} \right)
\sum_{n = 0}^{N - 1} e^{-j2 \pi f n} e^{j2 \pi n \frac{k}{N}} = \sum_{n = 0}^{N-1} e^{-j2 \pi n \left( f - \frac{k}{N} \right)} = \frac{1 - e^{-j2 \pi N \left( f - \frac{k}{N} \right)}}{1 - e^{-j2 \pi \left( f - \frac{k}{N} \right)}} = \frac{e^{-j \pi N \left( f - \frac{k}{N} \right)} \left( e^{j \pi N \left( f - \frac{k}{N} \right) } - e^{-j \pi N \left( f - \frac{k}{N} \right)} \right)}{e^{-j \pi  \left( f - \frac{k}{N} \right)} \left( e^{j \pi  \left( f - \frac{k}{N} \right) } - e^{-j \pi  \left( f - \frac{k}{N} \right)} \right)} = e^{-j \pi \left( N - 1 \right) \left( f - \frac{k}{N} \right)} \frac{\sin{\left( \pi N f - k \pi \right)}}{\sin{\left( \pi f - \pi \frac{k}{N} \right)}}

12-13 In pratica la DTFT viene approssimata a una DFT definita su un nuovo periodo N_1 \gg N, applicata a una sequenza x_z \left( n \right) che non è altro che la sequenza di partenza x \left( n \right) a cui sono stati accodati N_1 - N campioni nulli:

x_z \left( n \right) = \begin{cases} x \left( n \right) & n = 0 , \ldots , N - 1 \\
0 & n  = N , \ldots , N_1 - 1 \end{cases}

Aggiungere degli zeri in fondo alla sequenza x \left( n \right) permette di aumentare artificialmente il periodo della sua DFT senza alterare la DTFT a cui si cerca di approssimare:

X \left( e^{j2 \pi f_k} \right) = \sum_{n = 0}^{N-1} x \left( n \right) e^{-j2 \pi f_k n} = \sum_{n = 0}^{N_1 - 1} x_z \left( n \right) e^{-j 2 \pi \frac{k}{N_1} n} \quad k = 0, \ldots, N_1 - 1

e siccome il numero di campioni presi dalla DTFT è maggiore, la risoluzione della DTFT risulta molto più alta.

Esempio - Sequenza porta
x \left( n \right) = \begin{cases} 1 & 0 \leq n \leq N - 1 \\
0 & n \geq N \end{cases}

14 Considerando solamente l'intervallo \left[ 0 , N - 1 \right] si ottiene come DFT \left| X \left( k \right) \right| una funzione \mathrm{sinc} campionata a una frequenza molto bassa:

X \left( k \right) = \sum_{n = 0}^{N-1} x \left( n \right) e^{-j2 \pi n \frac{k}{N}} = \sum_{n = 0}^{N-1} e^{-j2 \pi n \frac{k}{N}} = \frac{1 - e^{-j2 \pi k}}{1- e^{-j2 \pi \frac{k}{N}}} = \begin{cases} 0 & k \neq 0 \\
N & k = 0 \end{cases}
N_1 = N
DFT zero-padding technique example (N) - Sequence.pdf   DFT zero-padding technique example (N) - DFT.pdf

16 Basta però aggiungere degli zeri a destra della porta per migliorare la risoluzione della DTFT:

N_1 = 2 N
DFT zero-padding technique example (2N) - Sequence.pdf   DFT zero-padding technique example (2N) - DFT.pdf

N_1 = 4 N
DFT zero-padding technique example (4N) - Sequence.pdf   DFT zero-padding technique example (4N) - DFT.pdf

19 N_1 = 100 N
DFT zero-padding technique example (100N) - DFT.pdf

Aliasing nel tempo[modifica]

22 Partendo dalla sequenza x \left( n \right) di durata non finita:

x \left( n \right) = a^n u \left( n \right) \quad 0 < a < 1

23 la sua DTFT campionata nel primo periodo N:

X \left( k \right) = X \left( e^{j2 \pi \frac{k}{N}} \right) = \frac{1}{1 - a^{-j2 \pi \frac{k}{N}}} \quad  k = 0, \ldots, N - 1

non coincide esattamente con la DFT calcolata sulla sequenza troncata:

X_N \left( k \right) = \frac{1 - a^N}{1 - a^{-j 2 \pi \frac{k}{N}}} \quad k = 0 , \ldots, N - 1

24 perché vi è un effetto di aliasing nel tempo dovuto al fatto che la sequenza x \left( n \right) di partenza non ha durata finita.

25 Al crescere dell'ampiezza del periodo N, la DFT della sequenza troncata tende sempre di più alla DTFT campionata nel primo periodo.

Proprietà della DFT[modifica]

Accorgimenti[modifica]

30 Le proprietà della DFT sono analoghe a quelle della DTFT, ma ci vogliono alcuni accorgimenti:

  • la DFT X \left( k \right) è periodica di periodo N → le operazioni sulla DFT X \left( k \right) devono condurre a una sequenza periodica di periodo N;
  • la IDFT x \left( n \right) è periodica di periodo N → le operazioni sulla IDFT x \left( n \right) devono condurre a una sequenza periodica di periodo N.

Operatore di modulo[modifica]

32 L'operatore di modulo restituisce sempre un numero compreso in \left[ 0 , N - 1 \right]:

{\left| k \right|}_N = k \text{ mod } N

Se k è un numero negativo, si somma N a k tante volte fino a ottenere un numero compreso in \left[ 0 , N - 1 \right]:

{\left| - 13 \right|}_{16} = 3

Operatore di ritardo circolare[modifica]

DFT circular shift.png

L'operatore di ritardo circolare restituisce una sequenza ritardata che è ancora compresa in \left[ 0 , N - 1 \right]:

x \left( {\left| n - N_0 \right|}_k \right)

in quanto i campioni che vanno al di là del primo periodo ricompaiono all'inizio del primo periodo stesso.

Modo operativo
  • si genera la sequenza periodicizzata \overline{x} \left( n \right);
  • si applica il ritardo di N_0 campioni a \overline{x} \left( n \right): \overline{x} \left( n - N_0 \right);
  • la sequenza ritardata x \left( n - N_0 \right) è pari alla sequenza periodicizzata \overline{x} \left( n  - N_0 \right) troncata nel primo periodo N (n \in \left[ 0 , N - 1 \right]).

Convoluzione circolare[modifica]

34 La convoluzione circolare tra due sequenze x \left( n \right) e y \left( n \right) è definita:

x \left( n \right) \otimes y \left( n \right) =  \sum_{k = 0}^{N-1} x \left( k \right) y \left( {\left| n - k \right|}_N \right)

dove N è la più grande tra le durate delle singole sequenze.

Proprietà commutativa
x \left( n \right) \otimes y \left( n \right) = y \left( n \right) \otimes x \left( n \right)

39 La convoluzione circolare si può rappresentare con un prodotto matrice per vettore. Ad esempio, se la durata N è pari a 4:

x \left( n \right) \otimes y \left( n \right) = \begin{bmatrix} y \left( 0 \right) & y \left( 3 \right) & y \left( 2 \right) & y \left( 1 \right) \\
y \left( 1 \right) & y \left( 0 \right) & y \left( 3 \right) & y \left( 2 \right) \\
y \left( 2 \right) & y \left( 1 \right) & y \left( 0 \right) & y \left( 3 \right) \\
y \left( 3 \right) & y \left( 2 \right) & y \left( 1 \right) & y \left( 0 \right)  \end{bmatrix} \begin{bmatrix} x \left( 0 \right) \\ x \left( 1 \right) \\ x \left( 2 \right) \\ x \left( 3 \right) \end{bmatrix}

Linearità[modifica]

z \left( n \right ) = a_1 x \left( n \right) + a_2 y \left( n \right) \Rightarrow Z \left( k \right) = a_1 X \left( k \right) + a_2 Y \left( k \right)

Ritardo[modifica]

z \left( n \right) = x \left( {\left| n - N_0 \right|}_N \right) \Rightarrow Z \left( k \right)  = X \left( k \right) e^{-j2 \pi \frac{k}{N} N_0}

Modulazione[modifica]

z \left( n \right) = e^{j2 \pi \frac{k_0}{N} n} x \left( n \right) \Rightarrow Z \left( k \right) = X \left( {\left| k - k_0 \right|}_N \right)

Convoluzione circolare[modifica]

33 z \left ( n \right)  = x \left( n \right) \otimes y \left( n \right) \Rightarrow Z \left( k \right) = X \left( k \right) \cdot Y \left( k \right)
Accorgimenti
  • la convoluzione circolare di due sequenze di durata N campioni ha una uguale durata di N campioni (a differenza della convoluzione lineare della DTFT, la cui sequenza risultante ha durata 2 N - 1 campioni);
  • la convoluzione circolare \overline{z} \left( k \right) di due IDFT \overline{x} \left( k \right) e \overline{y} \left( k \right), entrambe di periodo N, è periodica di periodo N:
    \overline{z} \left( n + N \right) = \sum_{k = 0}^{N-1} \overline{x} \left( k \right) \overline{y} \left( n + N - k \right) = \sum_{k = 0}^{N-1} \overline{x} \left( k \right) \overline{y} \left( n - k \right) = \overline{z} \left( n \right)

Proprietà di simmetria[modifica]

41 Se la sequenza x \left( n \right) è reale:

  • X \left( k \right) = X_R \left( k \right) + j X_I \left( k \right)
    • se x \left( n \right) è una sequenza pari: X_I \left( k \right) = 0;
    • se x \left(  n \right) è una sequenza dispari: X_R \left( k \right) = 0;
  • per la DFT vale la simmetria hermitiana intorno a 0:
    X \left( k \right) = X^* \left( {\left| - k \right|}_N \right)
  • 42 per la DFT vale la simmetria hermitiana intorno a \tfrac{N}{2};
  • il campione in k = 0 della DFT è sempre reale:
    X^* \left( 0 \right) = X \left( 0 \right) = \sum_{n = 0}^{N-1} x \left( n \right)

Trasformata di Fourier veloce (FFT)[modifica]

7-46 La DFT ha complessità quadratica (N^2).

47 La FFT è un algoritmo di complessità linearitmica (N \log{N}) che implementa la DFT e la IDFT in maniera più efficiente.

DFT[modifica]

50 La DFT può essere rappresentata in termini di matrici:

\mathbf{X} = \mathbf{H x} ; \; \begin{bmatrix} X \left( 0 \right) \\ X \left( 1 \right) \\ X \left( 2 \right) \\ \vdots \\ X \left( N - 1 \right) \end{bmatrix} = \begin{bmatrix} 1 & 1 & 1 & \cdots & 1 \\ 1 & H_N^1 & H_N^2 & \cdots & H_N^{N-1} \\ 1 & H_N^2 & H_N^4 & \cdots & H_N^{2 \left( N-1 \right)} \\
\vdots & \vdots & \vdots & \ddots & \vdots \\
1 & H_N^{N-1} & H_N^{2 \left( N - 1 \right)} & \cdots & H_N^{\left( N - 1 \right) \left( N - 1 \right)} \end{bmatrix} \cdot \begin{bmatrix} x \left( 0 \right) \\ x \left( 1 \right) \\ x \left( 2 \right) \\ \vdots \\ x \left( N - 1 \right) \end{bmatrix}

dove:

H_N = e^{-j 2 \pi \frac{1}{N}} \Rightarrow H_N^{nk} = e^{-j2 \pi n \frac{k}{N}}

Algoritmo della decimazione nel tempo[modifica]

52 L'algoritmo della decimazione nel tempo è un algoritmo di tipo "divide et impera" che riduce la complessità dell'algoritmo di DFT.

Ipotesi
N è una potenza di 2
Divide[modifica]

La sequenza si suddivide in due sottosequenze costituite da metà campioni:

  • la sottosequenza x_0 \left( n \right) è formata dai campioni pari:
    x_0  \left( n \right) = x \left( 2 n \right) \quad  n = 0 , \ldots , \frac{N}{2} - 1
  • la sottosequenza x_1 \left( n \right) è formata dai campioni dispari:
    x_1  \left( n \right) = x \left( 2 n + 1 \right) \quad  n = 0 , \ldots , \frac{N}{2} - 1
Ricombinazione[modifica]

53 Ogni DFT X \left( k \right), elemento del vettore \mathbf{x}, si ottiene dalla somma delle DFT delle singole sottosequenze x_0 \left( n \right) e x_1 \left( n \right):

X \left( k \right) = X_0 \left( {\left| k \right|}_{\frac{N}{2}} \right) + H_N^k X_1 \left( {\left| k \right|}_{\frac{N}{2}} \right) = \sum_{n = 0}^{\frac{N}{2} - 1} x_0 \left( n \right) H_{\frac{N}{2}}^{nk} + H_N^k \sum_{n = 0}^{\frac{N}{2} - 1} x_1 \left( n \right) H_{\frac{N}{2}}^{nk}
Dimostrazione
X \left( k \right) = \sum_{n = 0}^{N-1} x \left( n \right) H_N^{nk} = \sum_{n = 0}^{\frac{N}{2} - 1} x_0 \left(  n \right) H_N^{2nk} + \sum_{n = 0}^{\frac{N}{2} - 1} x_1 \left( n \right) H_N^{ \left( 2n+1 \right) k } = \sum_{n = 0}^{\frac{N}{2} - 1} x_0 \left(  n \right) H_N^{2nk} + H_N^k \sum_{n = 0}^{\frac{N}{2} - 1} x_1 \left( n \right) H_N^{  2n k } =

Siccome H_N^2 = H_{\frac{N}{2}}:

= \sum_{n = 0}^{\frac{N}{2} - 1} x_0 \left(  n \right) H_{\frac{N}{2}}^{nk} + H_N^k \sum_{n = 0}^{\frac{N}{2} - 1} x_1 \left( n \right) H_{\frac{N}{2}}^{  n k }
Complessità[modifica]

55 Il numero totale di operazioni (somme e prodotti) complesse è pari a:

N + {N^2 \over 2} < N^2

57 Si può quindi riapplicare ripetutamente il procedimento "divide et impera" sulle sottosequenze. Al passo k-esimo si ha una complessità pari a:

k \cdot N + 2^k {\left( \frac{N}{2^k} \right)}^2 \quad k = 1 , \ldots , \log_2 N

All'ultimo passo (k = \log_2 N), quando si arriva a dover valutare le DFT su 1 punto, la complessità totale è linearitmica:

\log_2 N \cdot N + N \approx \log_2 N \cdot N

Riduzione di complessità[modifica]

58 La complessità può essere ulteriormente ridotta sfruttando le proprietà di simmetria della matrice \mathbf{H}.

59-60 Esempio con N=4 campioni
\begin{cases} X \left( 0 \right) = x \left( 0 \right) + x \left( 2 \right) + x \left( 1 \right) + x \left( 3 \right) \\ X \left( 1 \right) =  x \left( 0 \right) + H_2 x \left( 2 \right) + H_4 \left[ x \left( 1 \right) + H_2 x \left( 3 \right) \right] \\ X \left( 2 \right) = x \left( 0 \right) + x \left( 2 \right) + H_4^2 \left[ x \left( 1 \right) + x \left( 3 \right) \right] \\ X \left( 3 \right) = x \left( 0 \right) + H_2 x \left( 2 \right) + H_4^3 \left[ x \left( 1 \right) + H_2  x \left( 3 \right) \right] \end{cases}
Dimostrazione
X \left( k \right) = X_0 \left( {\left| k \right|}_{\frac{N}{2}} \right) + H_N^k X_1 \left( {\left| k \right|}_{\frac{N}{2}} \right) = X_0 \left( {\left| k \right|}_2 \right) + H_4^k X_1 \left( {\left| k \right|}_2 \right) \quad k = 0, \ldots , 3
\begin{cases} X_0 \left( k \right) = x \left( 0 \right) + x \left( 2 \right) H_2^k \\ X_1 \left( k \right) = x \left( 1 \right)  + x \left( 3 \right) H_2^k \end{cases} \Rightarrow \begin{cases} X_0 \left( 0 \right) = x \left( 0 \right) + x \left( 2 \right) \\ X_0 \left( 1 \right) = x \left( 0 \right) + x \left( 2 \right) H_2 \\ X_1 \left( 0 \right) = x \left( 1 \right) + x \left( 3 \right) \\ X_1 \left( 1 \right) = x \left( 1 \right) + x \left( 3 \right) H_2 \end{cases} \Rightarrow \begin{cases} X \left( 0 \right) = X_0 \left( 0 \right) + X_1 \left( 0 \right) = x \left( 0 \right) + x \left( 2 \right) + x \left( 1 \right) + x \left( 3 \right) \\ X \left( 1 \right) = X_0 \left( 1 \right) + H_4 X_1 \left( 1 \right) = x \left( 0 \right) + H_2 x \left( 2 \right) + H_4 \left[ x \left( 1 \right) + H_2 x \left( 3 \right) \right] \\ X \left( 2 \right) = X_0 \left( 0 \right) + H_4^2 X_1 \left( 0 \right) = x \left( 0 \right) + x \left( 2 \right) + H_4^2 \left[ x \left( 1 \right) + x \left( 3 \right) \right] \\ X \left( 3 \right) = X_0 \left( 1 \right) + H_4^3 X_1 \left( 1 \right) = x \left( 0 \right) + H_2 x \left( 2 \right) + H_4^3 \left[ x \left( 1 \right) + H_2  x \left( 3 \right) \right] \end{cases}
61 Schema circuitale a farfalla

IDFT[modifica]

49 Tutte le considerazioni sulla DFT valgono anche per la IDFT: basta applicare l'algoritmo sul complesso coniugato della DFT e poi dividere per N:

x \left( n \right) = \frac{1}{N} \sum_{k = 0}^{N - 1} X \left( k \right) e^{j2 \pi n \frac{k}{N}} = \frac{1}{N} {\left[ \sum_{k=0}^{N-1} X^* \left( k \right) e^{-j2 \pi n \frac{k}{N}} \right]}^* \Rightarrow x \left( n \right) =  \text{IDFT} \left[ X \left( k \right) \right] = \frac{1}{N} \text{DFT} \left[ X^* \left( k \right) \right] , \quad  n = 0,1,2, \ldots, N - 1

Note[modifica]

  1. Si noti che la sequenza di partenza x \left( n \right) non era periodica, ma aveva durata finita.
Nota disambigua.svg Il titolo di questa pagina non è corretto per via delle caratteristiche del software MediaWiki. Il titolo corretto è B4. DFT e FFT.
Gtk-go-back-ltr.svg  B3. Trasformata di Fourier a tempo discreto Gtk-home.svg  Teoria ed elaborazione dei segnali (Elaborazione numerica dei segnali) B5. Analisi in frequenza di segnali a tempo continuo campionati  Gtk-go-forward-ltr.svg