Interpolazione trigonometrica in forma esponenziale

[dissonance]

Ho una funzione [tex]f\colon \mathbb{R}\to \mathbb{R}[/tex] periodica di periodo finito [tex]T[/tex]. Ne prendo [tex]N[/tex] campioni

[tex]$y_k=f(\frac{T}{N}k)[/tex] e ne faccio la DFT [tex]$Y_j=\frac{1}{N}\sum_{k=0}^{N-1}y_ke^{\imath(\frac{2\pi}{T}j)k}[/tex],

ottenendo il vettore [tex]\bold{Y}[/tex] di [tex]\mathbb{C}^N[/tex]. Ora come posso, tramite questo vettore, costruire il polinomio trigonometrico di interpolazione della [tex]f[/tex] nei nodi [tex]x_k=k\dfrac{T}{N}[/tex] in forma esponenziale (intendo una rappresentazione tipo [tex]$P(t)=\sum_{k=-a}^bc_ke^{\imath\frac{2\pi}{T}bt}[/tex])? Sul libro che sto leggendo si costruisce solo la rappresentazione con seni e coseni nel caso particolare [tex]T=2\pi[/tex] e mi sto incasinando con i conti.

[dissonance]

Ok, trovato. Che fatica, però! Riporto qui il risultato nel caso servisse a qualcuno. Per la definizione di [tex]f, \bold{y}, \bold{Y}, N, T[/tex] vedere il post precedente.

L'idea è quella di definire un polinomio trigonometrico assegnando ad ogni entrata del vettore [tex]\bold{Y}[/tex] una frequenza, secondo questo schema:

1)[tex]N[/tex] è dispari, [tex]N=2M-1[/tex], allora si suddivide [tex]\bold{Y}[/tex] come segue:

[tex]Y_0, \left[ Y_1 \ldots Y_{M-1}\right],\ \left[Y_M \ldots Y_{N-1} \right][/tex]

assegnando ai coefficienti del primo gruppo le armoniche corrispondenti, e ai coefficienti del secondo gruppo le armoniche negative secondo la corrispondenza: a [tex]Y_{N-j}[/tex] corrisponde la [tex]-j[/tex]-esima armonica.

2) [tex]N[/tex] è pari, [tex]N=2M[/tex], allora si suddivide [tex]\bold{Y}[/tex] come segue:

[tex]Y_0, \left[ Y_1 \ldots Y_{M-1}\right],\ Y_M,\ \left[Y_{M+1} \ldots Y_{N-1} \right][/tex]

poi si procede come nel caso 1).



***Formalmente***:



Sia [tex]\omega_N=e^{\imath \frac{2\pi}{N}}[/tex]. E' un fatto algebrico noto che [tex]\omega_N[/tex] è radice primitiva N-esima dell'unità; in particolare l'aritmetica delle sue potenze è modulo [tex]N[/tex]:

(*) [tex]\omega_N^k=\omega_N^h\quad \iff k \equiv h \mod N[/tex].

Per costruzione il polinomio [tex]$P(w)=\sum_{j=0}^{N-1}Y_j w^j[/tex] soddisfa le condizioni

[tex]P(\omega_N^k)=y_k,\quad k=0\ldots N-1[/tex];

da queste ricaviamo il polinomio trigonometrico di interpolazione della funzione [tex]f[/tex] distinguendo due casi:

1) [tex]N=2M-1[/tex]:
in questo caso possiamo riscrivere la condizione [tex]P(\omega_N^k)=y_k[/tex] come

[tex]$y_k=\sum_{j=0}^{N-1}Y_j\omega_N^{kj}=\sum_{j=0}^{M-1}Y_j\omega_N^{kj}+\sum_{j=M}^{2M-2}Y_j\omega_N^{kj}[/tex];

usando, nella seconda sommatoria, il cambiamento di indice [tex]h=N-j[/tex], e ricordando la (*), abbiamo allora

[tex]$y_k=\sum_{j=0}^{M-1}Y_j\omega_N^{kj}+\sum_{h=1}^{M-1}Y_{N-h}\omega_N^{k(-h)}[/tex].

Quindi il polinomio trigonometrico

[tex]$F(t)=\sum_{j=0}^{M-1}Y_je^{\imath(\frac{2\pi}{T}j)t}+\sum_{h=1}^{M-1}Y_{N-h}e^{\imath(\frac{2\pi}{T}(-h))t[/tex]

verifica le identità

[tex]$F(k\frac{T}{N})=y_k=f(k\frac{T}{N}),\quad k=0...N-1[/tex].

2) [tex]N=2M[/tex]
in questo caso considereremo a parte l'addendo corrispondente alla [tex]M[/tex]-esima frequenza:

[tex]$y_k=\sum_{j=0}^{M-1}Y_j\omega_N^{kj}+Y_M\omega_N^{kM}+\sum_{h=1}^{M-1}Y_{N-h}\omega_N^{k(-h)}[/tex].

Osserviamo poi che [tex]\omega_N^{kM}=\Re(\omega_N^{kM})=\cos(\pi k)[/tex].

Il polinomio trigonometrico cercato è

[tex]$F(t)=\sum_{j=0}^{M-1}Y_je^{\imath(\frac{2\pi}{T}j)t}+Y_M\cos(\frac{2\pi}{T}Mt)+\sum_{h=1}^{M-1}Y_{N-h}e^{\imath(\frac{2\pi}{T}(-h))t[/tex].

_________________________________
NOTA BENE: Anche se nel corso di questi post ho usato l'articolo determinativo (il polinomio trigonometrico di interpolazione), non è vero che esso è unico. La scelta indicata è però la migliore per vari motivi. Vedi http://en.wikipedia.org/wiki/Discrete_F ... polynomial

[dissonance]

Per concludere aggiungo la traduzione dei polinomi trigonometrici trovati sopra nella rappresentazione reale in seni e coseni. Le ipotesi sono quelle precedenti, in particolare ora diventa necessario che [tex]f[/tex] assuma valori reali.

Def.: Per ogni [tex]j=0\ldots N-1[/tex], siano [tex]$\frac{\alpha_j-\imath \beta_j}{2}=Y_j[/tex].
Risulta che

[tex]$\alpha_j= \frac{2}{N}\sum_{k=0}^{N-1}y_k\cos(\frac{2\pi}{N}jk)[/tex];

[tex]$\beta_j= \frac{2}{N}\sum_{k=0}^{N-1}y_k\sin(\frac{2\pi}{N}jk)[/tex].

________________

1)[tex]N=2M-1[/tex]

In questo caso è

[tex]$F(t)=Y_0+ \sum_{j=1}^{M-1}Y_je^{\imath (\frac{2\pi}{T}j)t}+\sum_{h=1}^{M-1}Y_{N-h}e^{\imath (\frac{2\pi}{T}(-h))t}[/tex];

ma osserviamo che [tex]Y_{N-h}=\overline{Y_h}[/tex] per ogni [tex]h=1\ldots M-1[/tex] perché il vettore [tex]\bold{y}[/tex] è reale; dunque

[tex]$F(t)=Y_0+\sum_{j=1}^{M-1}2\Re(Y_je^{\imath (\frac{2\pi}{T}j)t})[/tex]

e a conti fatti

[tex]$F(t)=\frac{\alpha_0}{2}+\sum_{j=1}^{M-1}\left( \alpha_j\cos(\frac{2\pi}{T}jt) + \beta_j \sin(\frac{2\pi}{T}jt) \right)[/tex].

2)[tex]N=2M[/tex]

In questo caso è

[tex]$F(t)=Y_0+ \sum_{j=1}^{M-1}Y_je^{\imath (\frac{2\pi}{T}j)t}+Y_M\cos(\frac{2\pi}{T}Mt)+\sum_{h=1}^{M-1}Y_{N-h}e^{\imath (\frac{2\pi}{T}(-h))t}[/tex];

si arriva, con ragionamento analogo a quello del caso 1), al risultato

[tex]$F(t)=\frac{\alpha_0}{2}+\sum_{j=1}^{M-1}\left( \cos(\frac{2\pi}{T}j t)+\beta_j \sin(\frac{2\pi}{T}j t) \right) + \frac{\alpha_M}{2}\cos(\frac{2\pi}{T}M t)[/tex].

L'unico punto da osservare esplicitamente è l'identità

[tex]$Y_M=\frac{\alpha_M}{2}[/tex]; questo è immediato se si osserva che, essendo

[tex]$Y_M=\frac{1}{N}\sum_{k=0}^{N-1}y_k e^{-\imath (\frac{2\pi}{N}\frac{N}{2})k}[/tex]

esso è un numero reale.