B.1 Caso finito
Fissano \(n \in \N\), si consideri un segnale definito (o misurato) su un intervallo discreto di \(n\) valori \[g : \cur{0, 1, \ldots, (n-1)} \to \mathbb{C}, \quad t \mapsto g( t ).\] Si definisce la sua trasformata di Fourier come la funzione \[ \hat g : \cur{0, 1, \ldots, (n-1)} \to \mathbb{C}, \quad \xi \mapsto \hat{g}(\xi) := \sum_{t=0}^{n-1} g(t) e^{ - 2 \pi i t \xi/n}.\]
Remark. Il dominio di definizione della \(g\) può essere pensato come un insieme di tempi, mentre il dominio di \(\hat{g}\) quello di opportune frequenze \(\xi\) per le funzioni oscillanti \(t \mapsto e^{ - 2 \pi i t \xi/n}\). Precisamente, le frequenze sarebbero \(\xi/n\), mentre le frequenze angolari \(2 \pi \xi/n\): questo giustifica parametrizzazioni diverse della trasformata di Fourier, ma nel caso discreto quella introdotta sopra è la più comune.
Se si interpeta sia \(g = (g(t))_{t=0}^{n-1}\) che \((\hat g(\xi))_{\xi =0}^{n-1}\) come vettori in \(\mathbb{C}^n\), allora \(\hat{g}\) è il vettore ottenuto moltiplicando \(g\) per la matrice \(F \in \mathbb{C}^{n\times n}\), data da \[ F_{\xi t } = e^{ - 2 \pi i t \xi/n}.\] La proprietà fondamentale della matrice \(F\) è di essere unitaria, ossia l’inversa di \(F\) è la sua trasposta coniugata (in realtà, per via della definizione che abbiamo usato, questo è vero a meno di una costante moltiplicativa \(1/n\)). Questo perché vale la relazione di ortogonalità, per \(s\), \(t \in \cur{0,1, \ldots, n-1}\), \[ (\bar{F}^T F)_{s t} = \sum_{\xi = 0}^{n-1} e^{ 2 \pi i s \xi/n}e^{ - 2 \pi i t \xi/n} = \begin{cases} n & \text{se $s=t$}\\ 0 &\text{altrimenti.} \end{cases}\] Per dimostrarlo basta ricordare la somma geometrica \(\sum_{j=0}^{k-1} z^j = (z^{k}-1)/(z-1)\) e il fatto che \(e^{2\pi i }= 1\).
Come prima conseguenza otteniamo allora la formula di inversione \[ g = \frac{1}{n} \bar{F}^T F g,\] che esplicitamente diventa \[ g(t ) = \frac 1 n \sum_{\xi=0}^{n-1} \hat{g}(\xi) e^{2 \pi i \xi t/n}.\] In altre parole, la trasformata di Fourier permette di ricostruire esattamente \(g\) mediante una operazione inversa che è analoga a quella diretta.
Una seconda conseguenza è il fatto che la norma (Euclidea) del vettore \(g\) coincide (a meno di un fattore \(1/n\)) con quella del vettore \(\hat{g}\), perché \[ |\hat g|^2 = \overline{F g}^T F g = \bar g^T \bar{F}^T F g = n \bar{g}^T g = n |g|^2.\]
Remark. La norma \(|g|^2\) può essere intepretata come una energia del segnale \(g\), di conseguenza l’identità sopra mostra che la stessa energia può essere ottenuta sommando le energie associate alle singole frequenze, ossia \(|\hat{g}(\xi)|^2\) (e dividendo per \(n\)).
Remark. Tutte le trasformate di Fourier che si approssimano numericamente sono ridotte al caso di tempi finiti. Per questo vi sono algoritmi particolarmente veloci, che in R si possono usare mediante la funzione fft()
. Ecco un esempio.
# fissiamo n
<- 16
n
# definiamo un vettore dei tempi e uno
# delle frequenze
<- 0:(n - 1)
t <- 0:(n - 1)
xi
# questa opzione permette di
# visualizzare due grafici uno accanto
# all'altro (1,2)=(1 riga, 2 colonne)
par(mfrow = c(1, 2))
# definiamo g come l'onda quadra e la
# plottiamo
<- c(rep(1, n/2), rep(-1, n/2))
g
plot(t, g, col = miei_colori[1], lwd = 3,
pch = 16)
# usiamo fft() per calcolare la
# trasformata di Fourier e la plottiamo
<- fft(g)
hat_g
<- max(abs(hat_g))
m
plot(xi, Re(hat_g), col = miei_colori[2],
ylab = "trasformata di g", pch = 16,
lwd = 3, ylim = c(-m, m))
points(xi, Im(hat_g) + 0.1, col = miei_colori[3],
pch = 16, lwd = 3)
legend("bottomright", c("parte reale", "parte immaginaria"),
fill = miei_colori[2:3], cex = 0.7)
Possiamo anche verificare la formula di inversione, usando l’opzione inverse =TRUE
nella stessa funzione fft()
. Bisogna tuttavia ricordare il fattore \(n\).
<- fft(hat_g, inverse = TRUE)
g_ricostruita
# osserviamo che la ricostruzione
# coincide con la g ma dilatata di un
# fattore n
plot(t, g_ricostruita, pch = 16, lwd = 3,
col = miei_colori[4])
Possiamo infine verificare l’identità dell’energia calcolando e confrontando le norme Euclidee dei vettori:
sum(abs(g)^2)
## [1] 16
sum(abs(hat_g)^2)
## [1] 256
# moltiplicando per n la prima si
# ottiene la seconda, infatti
* sum(abs(g)^2) n
## [1] 256