Interpolación trigonométrica usando la Transformada Discreta de Fourier

Introducción


La Transformada Discreta de Fourier (TDF) transforma una sucesión de $n$ números complejos 
$$\{x_0, x_1, \ldots, x_{N-1}\}$$
en otra sucesión de números complejos
$$\{X_0, X_1, \ldots, X_{N-1}\}$$
por medio de la fórmula
\begin{eqnarray} X_k &=&\; \frac{1}{N} \sum_{n=0}^{N-1} x_n \cdot e^{- i \, k \frac{2 \pi}{N}\, n}\qquad(k = 0, 1,\ldots, N-1)\nonumber\\ &=& \; \frac{1}{N} \sum_{n=0}^{N-1} x_j \left[ \cos \left( k \, \frac{2\pi}{N} \, n \right) - i\, \sin \left(k\, \frac{2\pi}{N}\,n\right) \right]\label{four} \end{eqnarray}
La segunda expresión se obtuve usando la fórmula de Euler $e^{it} = \cos t + i \text{sin}\, t$.

Veamos un ejemplo sencillo. Consideremos la siguiente sucesión de números complejos:
$$\mathbf x = \left\{ 1, 2-i, -i, -1+2i \right\}$$

Al aplicar la TDF (\ref{four}) tenemos:
\begin{equation} \begin{aligned} X_{0}= {} &\frac{1}{4} \big[ 1\cdot e^{-i\,2\pi \cdot 0\cdot 0/4}+\left(2-i\right) \cdot e^{-i\,2\pi \cdot 0\cdot 1/4}\big.\\ & \left.+\; \left(-i\right) \cdot e^{-i\,2\pi \cdot 0\cdot 2/4}+ \left(-1+2i\right) \cdot e^{-i\,2\pi \cdot 0\cdot 3/4}\right]=\frac{1}{4}(2)=\frac{1}{2}\\ X_{1}= {} & \frac{1}{4} \left[ 1\cdot e^{-i\,2\pi \cdot 1\cdot 0/4}+\left(2-i\right) \cdot e^{-i\,2\pi \cdot 1\cdot 1/4}\right.\\ & \left.+\;\left(-i\right)\cdot e^{-i\,2\pi \cdot 1\cdot 2/4}+ \left(-1+2i\right) \cdot e^{-i\,2\pi \cdot 1\cdot 3/4}\right]=\frac{1}{4}(-2-2i)= -\frac{1}{2}-\frac{1}{2}i\\ X_{2} = {} & \frac{1}{4} \left[ 1\cdot e^{-i\,2\pi \cdot 2\cdot 0/4}+\left(2-i\right) \cdot e^{-i\,2\pi \cdot 2\cdot 1/4}\right.\\ & \left.+\;\left(-i\right)\cdot e^{-i\,2\pi \cdot 2\cdot 2/4}+ \left(-1+2i\right) \cdot e^{-i\,2\pi \cdot 2\cdot 3/4}\right]=\frac{1}{4}(-2i) = -\frac{1}{2}i\\ X_{3} = {} & \frac{1}{4}\left[ 1\cdot e^{-i\,2\pi \cdot 3\cdot 0/4}+\left(2-i\right) \cdot e^{-i\,2\pi \cdot 3\cdot 1/4}\right.\\ & \left.+\;\left(-i\right)\cdot e^{-i\,2\pi \cdot 3\cdot 2/4}+ \left(-1+2i\right) \cdot e^{-i\,2\pi \cdot 3\cdot 3/4}\right]=\frac{1}{4}(4+4i)=1+i \end{aligned} \end{equation}

De esta forma, obtenemos la sucesión
$$\mathbf X = \left\{ \frac{1}{2},  -\frac{1}{2}-\frac{1}{2}i, -\frac{1}{2}i, 1+i \right\}.$$
Esto lo podemos escribir de forma abreviada como
$$\mathcal F (\mathbf x) =\mathbf X.$$

Muy bien, a partir de una sucesión de números complejos dada $\mathbf  x$ hemos obtenido otra sucesión de números complejos $\mathbf  X$. Pero, ¿esto de qué nos sirve? Desafortunadamente, el ejemplo anterior no da pistas suficientes para poder responder esta pregunta. Lo que necesitamos aquí es dar contexto a los datos, es decir, necesitamos mencionar de dónde proviene la sucesión $\mathbf  x$ para poder dar significado a la sucesión $\mathbf  X$ obtenida por medio de la TDF. Para este propósito, analicemos el problema de aproximar una función real definida en un intervalo $[a,b]$ usando lo que se conoce como interpolación trigonométrica.

Aproximación de funciones reales con funciones trigonométricas


Consideremos el siguiente problema: Sea $f: [0, 2 \pi]\rightarrow \mathbb R$ una función definida por
\begin{equation}\label{func} f(x) = 2x - \frac{x^2}{\pi}. \end{equation}



Dados $n$ valores de la función $f$, ¿es posible aproximar $f$ por medio de un polinomio trigonométrico definido en el mismo intervalo $[0,2\pi]$?

Por supuesto, la respuesta es afirmativa. De hecho, este problema se puede resolver usando diferentes métodos de interpolación. En particular, usaremos la TDF para calcular un polinomio trigonométrico cuyos valores coincidirán en su mayoría con los valores de $f$.

Nuestro objetivo aquí es aproximar la función $f$ usando funciones trigonométricas $ \text{sen}$ y $ \cos$. Para ello usaremos la interpolación trigonométrica en el plano complejo, conocida también como suma de Fourier, para determinar $p(x)$ tal que
\begin{equation}\label{approxim} f(x) \approx p(x) = c_0 + c_1 e^{i\,x }+ c_2 e^{2\, i\,x } + \cdots +c_{n-1} e^{(n-1)\,i\,x } = \sum_{k=0}^{N-1}c_k e^{i\,k\,x }. \end{equation}
Los valores $c_k$ se denominan coeficientes de Fourier. 

Ahora, dado que nuestra función (\ref{func}) es real, solo necesitamos usar la parte real de $p(x)$, es decir
\begin{equation}\label{approxre} f(x) \approx p(x) = c_0 + c_1 \cos x+ c_2 \cos 2x + \cdots +c_{n-1} \cos (n-1)x. \end{equation}
donde $c_{j}$ son constantes reales y el símbolo $\approx$ significa que la función $f(x)$ y la suma $p(x)$ coinciden en los puntos:
$$f(x_j) = p(x_j)\quad j=0, 1, 2, \ldots n-1.$$

Para simplificar aún más el problema dividiremos el intervalo $0\leq x \leq 2\pi$ en partes iguales, es decir, usaremos la partición
$$x_0= 0, x_1=\frac{2\,\pi}{N}, x_2=\frac{4\,\pi}{N},\cdots, x_j=\frac{2\,j\,\pi}{N},\cdots ,x_{N-1}=\frac{2(N-1)\pi}{N}.$$
De esta forma, si conocemos $N$ valores de la función $f$
$$\mathbf x = \{f_0(x_0), f_1(x_1), \ldots, f_{N-1}(x_{N-1})\},$$
entonces tenemos que calcular los coeficientes de Fourier
$$\mathbf X = \{c_0, c_1, \ldots, c_{N-1}\}.$$
Es aquí donde entra en acción la magia de la Transformada Discreta de Fourier. 

Por ejemplo, para $N=4$, tenemos como datos iniciales
\begin{equation}\label{senal}\mathbf x = \left\{f_0(0), f_1\left(\frac{\pi}{2}\right), f_2\left(\pi\right), f_3\left(\frac{3\pi}{2}\right) \right\} = \left\{0, \frac{3\pi}{4}, \pi, \frac{3\pi}{4} \right\}.\end{equation}



Como en este caso estamos considerando valores reales, usaremos la parte real de la TDF definida en (\ref{four}), es decir
\begin{equation} \begin{aligned} c_{0}= {} &\frac{1}{4}\left[ 0 \cdot \cos \left(2\pi/4 \cdot 0\cdot 0\right)+\left(3\pi/4\right) \cdot \cos\left(2\pi/4 \cdot 0\cdot 1\right) \right.\\ & \left. +(\pi)\cdot \cos\left(2\pi/4 \cdot 0\cdot 2\right)+(3\pi/4)\cdot \cos\left(2\pi/4 \cdot 0\cdot 3\right)\right]=\frac{5\pi}{8}\\ c_{1}= {} & \frac{1}{4} \left[ 0\cdot \cos (2\pi/4 \cdot 1\cdot 0)+(3\pi/4) \cdot \cos(2\pi/4 \cdot 1\cdot 1) \right.\\ & \left. + (\pi)\cdot \cos(2\pi/4 \cdot 1\cdot 2)+ (3\pi/4) \cdot \cos(2\pi/4\cdot 1\cdot 3)\right]=-\frac{\pi}{4}\\ c_{2}= {} & \frac{1}{4}\left[0 \cdot \cos (2\pi/4 \cdot 2\cdot 0)+(3\pi/4) \cdot \cos(2\pi/4\cdot 2 \cdot 1) \right.\\ & \left. +(\pi)\cdot \cos(2\pi/4 \cdot 2\cdot 2) +(3\pi/4) \cdot \cos(2\pi/4 \cdot 2\cdot 3)\right]=-\frac{\pi}{8}\\ c_{3}= {} & \frac{1}{4} \left[0 \cdot \cos (2\pi/4 \cdot 3\cdot 0)+(3\pi/4) \cdot \cos(2\pi/4\cdot 3\cdot 1) \right.\\ & \left. +(\pi)\cdot \cos(2\pi/4 \cdot3\cdot 2)+ (3\pi/4)\cdot\cos(2\pi/4\cdot 3\cdot 3)\right]=-\frac{\pi}{4} \end{aligned} \end{equation}
Por lo tanto, la interpolación trigonométrica está dada por el polinomio
\begin{equation}\label{firstapprx} p(x) = \frac{5\pi}{8} -\frac{\pi}{4} \cos x -\frac{\pi}{8} \cos 2x -\frac{\pi}{4} \cos 3x. \end{equation}
Esto se puede apreciar en la siguiente imagen.



Observación: El proceso de encontrar los coeficientes de Fourier a partir de una señal ($\ref{senal}$) se denomina Trasformada Discreta de Fourier. Mientras que el proceso inverso de reconstruir una señal a partir de sus coeficientes de Fourier usando la suma ($\ref{approxim}$) es conocido como la Inversa de la Transformada Discreta de Fourier.


Analicemos ahora el mismo problema de aproximar la función pero considerando más datos, por ejemplo para $N=16$. Debido a que el número de operaciones será muy grande, usaremos GeoGebra o JavaScript para realizar los cálculos que manualmente tomarían demasiado tiempo. En el siguiente applet puedes explorar la aproximación para diferentes valores de $N$.


Como habrás notado en el applet anterior mientras más grande sea el valor de $N$ mayor será la oscilación de la función $p(x)$

Las gráficas resultantes señalan una dificultad significativa con la TDF definida en (\ref{four}). Si bien los polinomios trigonométricos coinciden correctamente con los valores de la función dada, su comportamiento oscilatorio pronunciado los hace completamente inadecuados para la interpolación lejos de los puntos donde coinciden ambas funciones.

Sin embargo, esta dificultad se puede rectificar de la siguiente manera. El problema es que no estamos tomando atención a las frecuencias que están representadas en la suma de Fourier (\ref{approxim}). Mientras que la primera mitad de los sumandos en (\ref{approxim}) representan frecuencias relativamente bajas, la segunda mitad no, y puede ser reemplazada por una frecuencia más baja equivalente y, por lo tanto, con exponenciales menos oscilatorios. A saber, si $0< k \leq \frac{1}{2}N$ entonces 
$$e^{-i\,k\,x} \quad \text{y} \quad e^{i\,(N-k)\,x}$$
tienen los mismos valores, pero el primero tiene una frecuencia más baja comparado con el segundo. De esta forma, para nuestros propósitos, reemplazaremos la segunda mitad de los sumandos en (\ref{approxim}) por su alternativa versión con baja frecuencia. Si $N= 2m+1$ es impar, entonces definimos
\begin{equation}\label{approxodd} \begin{aligned} \hat{p}(x) &= c_{-m} e^{-i\,m\,x}+ \cdots + c_{-1} e^{-ix} + c_0 + c_1 e^{i\,x} + \cdots + c_{m} e^{i\,m\,x} \\&= \sum_{k=-m}^{m} c_k e^{i\,k\,x}  \end{aligned} \end{equation}
como la interpolación de baja frecuencia. Si, por otra parte $N=2m$ es par, entonces
\begin{equation}\label{approxeven} \begin{aligned} \hat{p}(x) &= c_{-m} e^{-i\,m\,x}+ \cdots + c_{-1} e^{-ix} + c_0 + c_1 e^{i\,x} + \cdots + c_{m} e^{i\,(m-1)\,x}\\ &= \sum_{k=-m}^{m-1} c_k e^{i\,k\,x}. \end{aligned} \end{equation}
En este caso, usaremos la Transformada Discreta de Fourier de bajas frecuencias para el caso par $N = 2m $ definida como
\begin{eqnarray} X_k &=& \frac{1}{N} {\large\sum_{ n =0}^{N-1} x_{n} \cdot e^{- i \,k \,\frac{2 \pi}{N} \, n}},\quad (k = -N/2, \ldots, N/2-1) \nonumber\\ &=& \frac{1}{N} {\large\sum_{ n=0}^{N-1} }\, x_{n} \left[ \cos \left( \frac{2\pi}{n}\,k\,n\right) - i\, \sin \left( \frac{2\pi}{n}\,k\,n\right) \right]\label{TDF2} \end{eqnarray}

Al aplicar la parte real de (\ref{TDF2}) a los valores 
$$\mathbf x = \left\{0, \frac{3\pi}{4}, \pi, \frac{3\pi}{4} \right\},$$
obtenemos los coeficientes
$$\mathbf X = \left\{-\frac{\pi}{8}, -\frac{\pi}{4}, \frac{5\pi}{8}, -\frac{\pi}{4} \right\}.$$
De esta manera, podemos reemplazar (\ref{firstapprx}) con el polinomio trigonométrico
$$\hat{p}(x) = -\frac{\pi}{8}\cos\left(-2\cdot x\right)-\frac{\pi}{4}\cos\left(-1 \cdot x\right)+\frac{5\pi}{8}\cos\left(0\cdot x\right)-\frac{\pi}{4}\cos\left(x\right)$$
En el siguiente applet podemos apreciar que se trata de una mejor aproximación de la función $f$ para valores pares de $n$.


Lo mismo podemos hacer para valores impares de $n$. Usa el siguiente applet para explorar distintos valores.


Listo, hemos aproximado una función real definida en un intervalo usando la TDF (en su versión con bajas frecuencias). Lo más sorprendente aquí es que en realidad no necesitamos conocer la función $f$ desde el principio. De hecho, solo necesitamos conocer una cantidad finita de datos
$$\mathbf x = \{f_0, f_1, \ldots, f_{N-1}\},$$
para reconstruir dicha función. Claro que mientras mayor sea el número de datos, la aproximación será más precisa.


Además, no es necesario restringirnos a aproximar funciones en un solo intervalo $[a,b]$, en realidad podemos extender la aproximación a toda la recta real. En particular, la TDF es una herramienta excelente para aproximar funciones periódicas $f(x, t)$ (también llamadas funciones de onda). Considera una función $f(x,t)$ cuyo parámetro $t$ varía de forma constante como se muestra en la imagen siguiente:



Con la TDF podemos aproximar esta función usando funciones sinusoidales más simples, es decir, funciones $\text{sen}$ o $\cos.$ En otras palabras, la TDF descompone a una función dada $f$ en términos de funciones sinusoidales más simples. Esto se puede apreciar en el siguiente applet, usa el ratón para bajar el deslizador:


En efecto, $f(x,t)$ está definida como:
$$\text{sen}(x+t) + \text{sen} (3x+t)$$
En el applet anterior, la curva de en medio representa la función $\text{sen}(x+t) $ y curva de abajo representa la función $ \text{sen} (3x+t)$. 

Lo mismo podemos hacer para funciones de onda más interesantes. Consideremos, por ejemplo, la onda cuadrada.
Aunque no lo parezca, esta onda también se puede descomponer en funciones sinusoidales.
Sin embargo, en esta ocasión necesitamos muchos términos, técnicamente una cantidad infinita de sumandos para representarla perfectamente. A medida que sumamos más y más ondas sinusoidales, el patrón se acerca cada vez más a la onda cuadrada con la que comenzamos.

Visualmente, notarás que en realidad las primeras ondas sinusoidales son las que marcan la mayor diferencia. Con el control deslizante a la mitad, tenemos la forma general de la onda, pero todo es ondulante. Solo necesitamos el resto de los pequeños para hacer que la ondulación se aplane.

Conclusión

Con los ejemplos aquí mostrados hemos visto que la Transformada Discreta de Fourier nos permite transformar un conjunto finito de datos para aproximar y descomponer funciones periódicas en términos de funciones sinusoidales.

Esto es de gran utilidad en la actualidad debido al desarrollo de medios digitales y a la necesidad por transmitir información de forma eficiente. En los medios digitales modernos (audio, imágenes fijas o video), las señales continuas se muestrean a intervalos de tiempo discretos antes de ser procesadas. La Transformada Discreta de Fourier descompone señales muestra (que se pueden definir como funciones reales o complejas) en sus componentes periódicos fundamentales: senos y cosenos, o, más convenientemente, exponenciales complejos. El hecho crucial, en el que se basa todo el procesamiento moderno de señales, es que los exponenciales complejos muestreados forman una base ortogonal (tema que quizá exploraremos en otra entrada).

Referencias


Nota: Los applets en JavaScript fueron escritos originalmente por Jez Swanson. Todos los applets hechos con GeoGebra se pueden consultar en el libro en línea Interpolation using the Discrete Fourier Transform.

También te recomiendo la siguiente entrada para conocer otras divertidas aplicaciones de la Transformada Discreta de Fourier:




Gracias por llegar al final de este artículo. Si deseas puedes apoyarme en Patreon usando el siguiente enlace:

Become a Patron!

Con tu apoyo podré seguir escribiendo y compartiendo artículos y applets de matemáticas.

Comentarios

Entradas populares

Galileo Galilei y su ley de caída libre

Breve historia del Cálculo

Una historia de la Teoría de Conjuntos