```python
import holoviews as hv
hv.extension('bokeh')
```


```python
import numpy as np
import scipy.signal
```

# Espectro y Transformada de Fourier

## Espectro de frecuencia

### Breve historia

- Isaac Newton llamó **espectros** a los componentes que forman la luz blanca y que normalmente no se pueden ver
- Newton mostró usando prismas que la luz blanca puede descomponerse en colores y viceverza
- Hoy entendemos que la luz como onda tiene una **frecuencia asociada** y que cada color es una frecuencia particular


<table><tr><td>

  

</td><td>

  

</td></tr></table>


Paradójicamente, Newton nunca acepto que esto se debía a la frecuencia de la radiación ya que creía en la teoría corpuscular de la luz


### Ondas y señales sinusoidales


Una onda es una perturbación que transporta energía a través del espacio

Tipicamente se describe por su frecuencia, amplitud y desfase. La frecuencia (Hertz) es el recíproco del período (segundos)

La siguiente es una función del tiempo completamente descrita por su amplitud $A$, frecuencia $f$ (período $P=1/f$) y fase $\phi$

$$
s(t) = A \cos (2 \pi f t + \phi)
$$

Observe en el siguiente ejemplo como cambia la señal al modificar sus parámetros


```python
def wave(time, frequency, phase, amplitude):
    return amplitude*np.cos(2.0*np.pi*time*frequency + phase)
```


```python
signal_plot = hv.HoloMap(kdims=['Amplitud', 'Frecuencia', 'Desfase'])

t = np.linspace(0, 2, num=1000)
for A in [0.5, 1.]:
    for f in [0.5, 1, 2]:
        for phi in np.arange(0., np.pi, np.pi/10):
            signal_plot[(A, f, phi)] = hv.Curve((t, wave(t, f, phi, A)))
            
signal_plot.opts(hv.opts.Curve(width=500, height=300))
```

### Componentes frecuenciales y descomposiciones armónicas

Una sinusoide que es períodica en una fracción

$$
\frac{P}{k} ~ \forall k \in \mathbb{N}
$$ 

también es periódica en $P$

En general llamamos a 

- $f_0 = 1/P$ la **frecuencia fundamental**
- $f_k = k/P = kf_0$ el **k-esimo armónico de $f_0$**

Si sumamos armónicos con distintas amplitudes el resultado es una nueva función periódica que tiene la misma frecuencia fundamental

En el siguiente ejemplo modifique las amplitudes de cada armónico y observe la señal resultante

$$
s(t) = A_1 \cos(2\pi f_0 t ) + A_2 \cos(2\pi 2 f_0 t) + A_3 \cos(2\pi 3 f_0 t)
$$

Note que la suma de tres sinusoides no es necesariamente otra sinusoide


```python
t = np.linspace(0, 2, num=1000)

def my_signal(time, frequency, Amplitudes):
    s = 0
    for k, A in enumerate(Amplitudes):
        s += A*np.cos(2.0*np.pi*time*(k+1)*frequency) 
    return s
```


```python
signal_plot = hv.HoloMap(kdims=['Amplitud fundamental', 
                                'Amplitud 1er armónico', 
                                'Amplitud 2do armónico', 
                                'Frecuencia'])
for f in [1, 2]:
    for A1 in [0.5, 1.]:
        for A2 in [0, 0.5, 1.]:
            for A3 in [0, 0.5, 1.]:
                signal_plot[(A1, A2, A3, f)] = hv.Curve((t, my_signal(t, f, [A1, A2, A3])))
                
signal_plot.opts(hv.opts.Curve(width=500, height=300))
```

La noción de sintetizar señales en base a sinusoides se generaliza por medio de la **serie trigonométrica**

$$
\begin{align}
s(t) &= \sum_{k=0}^\infty A_k \cos(2\pi k f t + \phi_k) \nonumber \\
&= \sum_{k=0}^\infty a_k \cos(2\pi k f t) + b_k \sin(2\pi k f t), \nonumber
\end{align}
$$

donde $a_k = A_k \cos(\phi_k)$ y $b_k = -A_k \sin(\phi_k)$ se obtienen de $\cos(x+y) = \cos(x)\cos(y) - \sin(x)\sin(y)$

Con esto podemos generar funciones periódicas arbitrarias definiendo $\{a_k, b_k\}$ y $f$

**Ejemplo**

Sea por ejemplo $a_k = 0$ y $b_k = 1/k$  $\forall k$

¿Qué señal se obtiene al agregar cada vez más armónicos?


```python
t = np.linspace(0, 5, num=1000)

def trigonometric_series(time, frequency, K):
    s = 0
    for k in range(1, K):
        s += np.sin(2.0*np.pi*time*k*frequency)/(k) 
    return s
```


```python
signal_plot = hv.HoloMap(kdims=['Cantidad de armónicos', 'Frecuencia'])
for f in [1, 2]:
    for K in [1, 2, 3, 5, 10, 20, 50, 100]:     
        signal_plot[(K, f)] = hv.Curve((t, trigonometric_series(t, f, K)))
                
signal_plot.opts(hv.opts.Curve(width=500, height=300))
```

```{note}
Hasta ahora hemos visto como sintetizar señales a partir de sus armónicos, pero también podemos hacer el proceso inverso, es decir encontrar los armónicos de una señal dada
```    


Tal como Newton descompuso la luz blanca en colores, nosotros podemos descomponer una señal en sus armónicos usando la **serie de Fourier**

## Serie de Fourier

Un poco de historia: En 1807 *Jean Baptiste Joseph Fourier* presenta un teorema indicando que una función periódica arbitraria con periódo $P=1/f_0$ puede representarse como una suma ponderada de senos y cosenos


La serie de Fourier (FS) es una generalización de la serie trigonométrica a los números complejos. Puede revisar el apéndice de esta lección si necesita refrescar la memoria sobre los números complejos y sus representaciones

La FS se define como

$$
s(t) = \sum_{k=-\infty}^{\infty} c_k e^{j 2\pi k f_0 t}, ~~ c_k \in \mathbb{C}
$$


A continuación veremos como obtener los coeficientes de Fourier $c_k$ para una señal periódicas arbitraria

### Base de Fourier

El conjunto de funciones 

$$
v_k (t) = \frac{1}{\sqrt{P}} e^{j2\pi k t / P} ~~ \forall k \in \mathbb{Z}
$$

se conoce como base de Fourier y cumple con

$$
\langle v_n (t), v_m (t) \rangle = \int_0^P v_n (t) v_m^* (t) dt = \frac{1}{P} \int_0^P e^{j2\pi (n-m)t/P} dt =\begin{cases}1 & n=m \\ 0 & n \neq m\end{cases}
$$

La base de Fourier es un conjunto ortonormal en el espacio de funciones periódicas con periódo $P$. Esto es facilmente comprobable si estudiamos

$$
\int_0^P e^{j2\pi k t / P} dt = \int_0^P \cos(2\pi k t/P) dt + j \int_0^P \sin(2\pi k t/P) dt 
$$

Podemos usar esta propiedad y escribir el producto punto entre una señal y el m-esimo elemento de la base de Fourier

$$
\begin{align}
\langle s(t),  e^{j 2\pi m t/ P} \rangle &= \int_0^P s(t) e^{-j 2\pi m t/ P} dt \nonumber \\
&= \int_0^P \sum_{k=-\infty}^{\infty} c_k e^{j 2\pi k t/P} e^{-j 2\pi m t/P} dt \nonumber \\
&= \sum_{k=-\infty}^{\infty} c_k  \int_0^P e^{j 2\pi (k-m) t/P}  dt \nonumber \\
&=  c_m  P \nonumber \\
\end{align}
$$

Esto nos da una forma sencilla para obtener los $c_m$ de una señal periódica arbitraria

$$
c_m = \frac{1}{P} \int_0^P s(t) e^{-j 2\pi m t/P} dt
$$


### Ejemplo formativo: FS de señal cuadrada

Sea 

$$
s(t) = \begin{cases} 1 & t \in[0, \frac{P}{2}] \\ 0 & t \in [\frac{P}{2}, P] \end{cases}
$$

Los coeficientes de su FS son

$$
c_0 = \frac{1}{P} \int_0^P s(t) dt = \frac{1}{P} \int_0^{P/2} dt = \frac{1}{2}
$$

y

$$
\begin{align}
c_k &= \frac{1}{P} \int_0^\frac{P}{2} e^{-j2\pi kt/P} dt  \nonumber  \\  
&= -\frac{j}{P}  \int_0^\frac{P}{2} \sin(2\pi kt/P) dt \nonumber \\
&=  0 + j \frac{\cos(\pi k)  - 1}{2\pi k} \nonumber 
\end{align}
$$

Notemos que a excepción de $c_0$ los coeficientes sólo tienen parte imaginaria. Además sólo los armónicos impares (senos) son distintos de cero

Finalmente la FS está dada por

$$
\begin{align}
s(t) &= \sum_{k=-\infty}^{\infty} j \frac{\cos(\pi k)  - 1}{2\pi k}  e^{j 2\pi k t/P} \nonumber \\
&= \frac{1}{2} + \sum_{k=1}^{\infty}  \frac{1 - \cos(\pi k)}{\pi k} \sin(2\pi k t/P)  \nonumber 
\end{align}
$$

Programemos esta expresión y visualicemos el resultado al incrementar cada vez la cantidad de armónicos


```python
t = np.linspace(0, 5, num=1000)

def FS_square(tiempo, frecuencia, K):
    s = 0.5 # C_0
    for k in range(1, K):
        c_k = (1-np.cos(np.pi*k))/(np.pi*k)
        s += c_k*np.sin(2.0*np.pi*tiempo*k*frecuencia)
    return s
```


```python
signal_plot = hv.HoloMap(kdims=['Cantidad de armónicos', 'Frecuencia', ])

for f in [1, 2]:
    for K in [1, 2, 3, 5, 10, 20, 50, 100]:     
        signal_plot[(K, f)] = hv.Curve((t, FS_square(t, f, K)))
                
signal_plot.opts(hv.opts.Curve(width=500, height=300))
```

### Propiedades de la FS

A continuación se enumeran algunas de las propiedades más importantes de la FS

- Si $s(t)$ es par entonces $c_k$ es par
- Si $s(t)$ es impar entonces $c_k$ es impar
- Si $s(t + P/2) = -s(t)$ (antiperiódica) entonces $c_k=0$ para k par
- Si $s(t)$ es real y par entonces $c_k$ es real y par
- Si $s(t)$ es real e impar entonces $c_k$ es imaginario e impar
- **La FS es lineal**

El **Teorema de Parseval** relaciona la potencia de una señal periódica arbitraria con sus coeficientes de Fourier

$$
\begin{align}
P_s &= \frac{1}{P} \int_0^P |s(t)|^2 dt \nonumber \\
&= \frac{1}{P} \int_0^P |\sum_{k=-\infty}^{\infty} c_k e^{j 2\pi k t/P }|^2 dt \nonumber \\
&= \frac{1}{P} \int_0^P \sum_{k=-\infty}^{\infty} |c_k |^2 dt \nonumber \\
&= \sum_{k=-\infty}^\infty |c_k|^2 \nonumber
\end{align}
$$

## Transformada de Fourier

Podemos extender el concepto de descomposición armónica a señales no periódicas usando la **Transformada de Fourier**. Esta herramienta matemática será fundamental en este curso
    

- El concepto de frecuencia puede aplicarse también a señales no-periódicas
- Según **Joseph Fourier** una señal no-periódica puede ser vista como una señal periódica **con un período infinito**
- El único requisito es que ahora las frecuencias son un continuo, con un espaciado infinitesimal


A continuación veremos que una señal analógica puede ser vista como continua en el tiempo o continua en frecuencia. En la próxima lección veremos la Transformada de Fourier discreta (DFT) para señales digitales


### Derivación de la transformada de Fourier a partir de la FS

Sea un tren de pulsos cuadrado con periódo P y ancho $2T < P$ definido en un período como

$$
s(t) = \begin{cases} 1, & |t| < T \\ 0, & T<|t| < P/2 \end{cases}
$$

Los coeficientes de su serie de Fourier son

$$
c_0 = \frac{1}{P} \int_{-P/2}^{P/2} s(t) dt = \frac{1}{P} \int_{-T}^{T}  dt = \frac{2T}{P}
$$

y

$$
\begin{align}
c_k &= \frac{1}{P} \int_{-T}^{T} e^{-j2\pi kt/P} dt  \nonumber \\
&= \frac{1}{\pi k} \sin \left (2\pi k \frac{T}{P} \right)  \nonumber 
\end{align}
$$


En el siguiente ejemplo la figura izquierda muestra la señal de pulso cuadrado mientras que la derecha muestra $c_k$ en función de $k$

¿Qué ocurre con $c_k$ a medida que $P$ aumenta ($f$ disminuye)?


```python
def fs_coef(P, T, maxK=200):
    K = np.arange(-maxK, maxK)
    ck = np.zeros_like(K, dtype='float32')
    ck[len(K)//2] = 2*T/P
    for k in range(1, K[-1]):
        ck[len(K)//2+k] = np.sin(2*np.pi*k*T/P)/(np.pi*k)
        ck[len(K)//2-k] = ck[len(K)//2+k]
    return K, ck

def synthesis(time, ck, P):
    maxK = len(ck)//2
    s = ck[maxK]
    for k in range(1, maxK):    
        s += 2*np.cos(2*np.pi*time*k/P)*ck[maxK+k]
    return s   
```


```python
t = np.linspace(-3, 3, num=2000)
T = 0.25

signal_plot = hv.HoloMap(kdims=['Periodo'])
fs_plot = hv.HoloMap(kdims=['Periodo'])
for P in [1, 2, 3, 5, 10]:
    K, ck = fs_coef(P, T)
    fs_plot[P] =  hv.Spikes((K, P*ck), 'Frecuencia (k)', 'P x ck').opts(xlim=(-50, 50))
    signal_plot[P] =  hv.Curve((t, synthesis(t, ck, P)), 'Tiempo (segundos)', 'Señal')
    

(signal_plot+fs_plot)
```

De la figura tenemos que cuando P es grande (o equivalentemente f es pequeño) 

- el tren de pulsos tiende a un único pulso
- los coeficientes toman una "forma suave"


Llamaremos a esta forma suave de $P \cdot c_k$ la "envolvente" de los coeficientes. Si tomamos el caso $P \to \infty$ la envolvente resultante es

$$
\lim_{P\to \infty} P c_k = S(f) = \int_{-\infty}^{\infty} s(t) e^{-j 2\pi t  f} dt,
$$

o también

$$
S(\omega) = \int_{-\infty}^{\infty} s(t) e^{-j\omega t } dt,
$$

donde $\omega = 2\pi f$ se llama frecuencia angular.

```{note}
Esto se conoce como transformada de Fourier directa o integral de Fourier
```


Reemplazando $c_k = f_0 S(k f_0)$ y tomando el límite cuando $f_0 = \frac{1}{P} \to 0$ tenemos

$$
\begin{align}
s(t) &= \lim_{f_0 \to 0} \sum_{k=-\infty}^{\infty} f_0 S(k f_0) e^{j 2\pi t  k f_0} \nonumber  \\
&= \int_{-\infty}^{\infty} S(f) e^{j 2\pi t  f} df \nonumber 
\end{align}
$$

o también

$$
s(t) = \frac{1}{2\pi} \int_{-\infty}^{\infty} S(\omega) e^{j \omega t } d\omega, 
$$

que se conoce como la transformada de Fourier inversa




### Par de Fourier

Los operadores (transformada de Fourier directa) 

$$
S(\omega) = \mathbb{FT} [s(t)] = \int_{-\infty}^{\infty} s(t) e^{-j\omega t } dt
$$

y (transformada de Fourier inversa)

$$
s(t) = \mathbb{FT}^{-1} [S(\omega)] = \frac{1}{2\pi} \int_{-\infty}^{\infty} S(\omega) e^{j \omega t } d\omega, \nonumber 
$$

se conocen como par de Fourier y nos permiten analizar una señal en el dominio del tiempo o en el dominio de la frecuencia sin pérdidas

### Ejemplo: Transformada del pulso cuadrado

Consideremos nuevamente el pulso cuadrado

$$
s(t) = \begin{cases} 1, & |t| < T \\ 0, & |t| > T\end{cases}
$$

su transformada de Fourier es

$$
\begin{align}
S(\omega) &= \int_{-\infty}^{\infty} s(t) e^{-j\omega t } dt \nonumber \\
&= \int_{-T}^{T} e^{-j\omega t } dt \nonumber \\
&= \frac{1}{-j\omega} \left(e^{-j\omega T }  - e^{j\omega T } \right) \nonumber \\
&= \frac{2}{\omega}  \sin(\omega T) = 2T \text{sinc}(\omega T) \nonumber 
\end{align}
$$

Esto calza precisamente con la "envolvente" que vimos anteriomente


```python
T = 0.25

w = np.arange(-100, 100)
S = 2*np.sin(w*T)/w
S[w==0] = 0.5
```


```python
hv.Curve((w, S), 'Frecuencia angular', 'Espectro').opts(width=500)
```

### Propiedades de la transformada de Fourier


1) La transformada de Fourier es un operador lineal, si tenemos dos señales y dos valores escalares entonces

$$
\mathbb{FT}[c_1 s_1(t) + c_2 s_2(t)] = c_1\mathbb{FT}[s_1(t)] + c_2\mathbb{FT}[s_2(t)] 
$$

2) El operador de convolución en el tiempo se convierte en multiplicación en frecuencia 

$$
\mathbb{FT}[(s_1 * s_2)(t)] =  \mathbb{FT}[s_1(t)] \cdot \mathbb{FT}[s_2(t)],
$$

donde 

$$
(s_1 * s_2)(t) = \int s_1(\tau) s_2(t-\tau) d\tau
$$

es la operación de convolución

3) Así mismo, la multiplicación en frecuencia se convierte en multiplicación en el tiempo

$$
\mathbb{FT}[s_1(t)\cdot s_2(t)] =  \frac{1}{2\pi}\mathbb{FT}[s_1(t)] * \mathbb{FT}[s_2(t)]  
$$

Finalmente el Teorema de Parseval nos dice que la energía se preserva entre frecuencia y tiempo

$$
\int | s(t) |^2 dt = \frac{1}{2\pi} \int | S(\omega) |^2 d\omega
$$


### Definición: Espectros de amplitud y fase

Llamamos a $S(\omega)$ el espectro (transformada de Fourier) de una señal $s(t)$ 

El espectro es un número complejo que podemos escribir en notación polar como

$$
S(\omega) = |S(\omega)| e^{j\Phi(\omega)},
$$

donde $|S(\omega)|$ se conoce como **espectro de amplitud** y $\Phi(\omega)$ como **espectro de fase**

## Resumen de la lección

En esta lección aprendimos

- a realizar síntesis de señales usando sinusoides
- a descomponer señales periódicas usando la serie de Fourier
- a descomponer señales periódicas y no periódicas usando la transformada de Fourier



## Apéndice: Números complejos

Sea z un número complejo, lo podemos escribir en forma cartesiana 

$$
z = \Re[z] + j \Im[z] = a + j b
$$

donde $a \in \mathbb{R}$, $b \in \mathbb{R}$ y $j = \sqrt{-1}$ es el número imaginario.


También podemos escribirlo en forma polar

$$
z = c e^{j\phi} = c \cos(\phi) + j c \sin(\phi)
$$

donde 
- $c = |z| = \sqrt{a^2 + b^2} \in [0, \infty]$ es la magnitud 
- $\phi = \angle z = \tan^{-1} \left (\frac{b}{a} \right) \in [-\frac{\pi}{2}, \frac{\pi}{2}]$ es el ángulo
- $a = c \cos(\phi)$
- $b = c\sin(\phi)$

se pueden escribir las siguientes relaciones 

$$
\cos(\phi) = \frac{1}{2} (e^{j\phi} + e^{-j\phi}) ~\wedge~ \sin(\phi) = \frac{1}{2j} (e^{j\phi} - e^{-j\phi})
$$

el complejo conjugado de $z = a + j b = c e^{j\phi}$ es

$$
z^* = a - jb = c e^{-j\phi}
$$


```python

```

5a2929dc36a6f6d94b5ad6e65012899a87a94c14
ipynb
Jupyter Notebook
lectures/unit1/lecture2.ipynb
phuijse/UACH-INFO183
0e1b6bef0bd80cda2753bd11e62016268f2de638
MIT
2018-08-27T23:53:15.000Z
2021-08-16T23:31:05.000Z
lectures/unit1/lecture2.ipynb
phuijse/UACH-INFO183
0e1b6bef0bd80cda2753bd11e62016268f2de638
MIT
lectures/unit1/lecture2.ipynb
phuijse/UACH-INFO183
0e1b6bef0bd80cda2753bd11e62016268f2de638
MIT
2019-01-04T17:43:55.000Z
2021-09-07T16:07:18.000Z
Qwen/Qwen-72B
1. YES
2. YES