Usar la PC para ver espectros de señales
Código auxiliar
# Cargo paquetitos
using DSP, FFTW
# Y armo un par de funciones auxiliares
stem(args...; kwargs...) = sticks(args...;
marker=:circle,
leg=false,
kwargs...)
# Delta
d(n) = n == 0 ? 1. : 0.
# Escalón
u(n) = n >= 0 ? 1. : 0.
"""
Función módulo pero con offset (opcional)
Manda a `t` al intervalo [from, from+length)
sumándole o restándole algún múltiplo de `len`
"""
cshift(t, len; from=0) = mod(t - from, len) + from
using Plots
Plots.default(:legend, false)
Señales discretas
inline_formula not implemented para inline_formula not implemented
x1(n) = cos(2*π*n/10) * (u(n) - u(n-100));
La graficamos
let
ns = -5:105
xs = x1.(ns)
stem(ns, xs)
end
Queremos el espectro de esto, aunque sea aproximado, sin cuentas analíticas.
n_tot = 500 # la cantidad de muestras que voy a tomar. ¡Incluye todo el soporte!
ns = range(0; length=n_tot) # un vector de índices `n`
x1s = x1.(ns) # un vector con las primeras `n_tot` muestras
# Calculo su DFT
x1s_dft = fft(x1s)
# Armo un vector de las frecuencias donde caen las muestras
fs = range(0; step=2π/n_tot, length=n_tot)
# Grafico el módulo, uniendo las muestras con una recta
plot(fs, abs.(x1s_dft))
# Graficamos pero en el intervalo central
let
# Mandamos la segunda mitad al inicio (fftshift)
# ...y reacomodamos las frecuencias en el intervalo [-pi pi)
fs_aux = cshift.(fftshift(fs), 2π, from=-π)
plot(fs_aux, fftshift(abs.(x1s_dft)) )
end
Señal continua
Esto es lo mismo que antes, luego de un muestreo. Imaginemos una señal triangular que arranca en cero
x2(t) = tr(t/2 - 0.5)
# tr(t) es una funcióen triángulo centrado entre -0.5 y 0.5
function tr(t)
res = 1 - 2 * abs(t)
return -0.5 < t < 0.5 ? res : zero(res)
end
;
plot(x2, -1, 5)
La muestreamos con alguna frecuencia de muestreo `sr` lo suficientemente grande como para que la energía de la señal esté (casi) toda contenida en inline_formula not implemented.
sr = 100
x2d(n) = x2(n / sr)
;
y ahora es un problema como el de antes
Tomamos un número de muestras que contenga todo el soporte - i.e. mayor a inline_formula not implemented
DFT, y el resultado son muestras del espectro de la discreta
Lo grafico centrado
Ahora, como muestreamos "con Nyquist", podemos interpretar a las frecuencias discretas como frecuencias continuas inline_formula not implemented
Hagamos todo directo, sin pasar por la señal discreta.
# Los tiempos donde voy a muestrear para armar el vector
ts = range(0; stop = 3, step=1/sr)
x2s = x2.(ts)
x2s_fft = fft(x2s)
plot(
range(-sr/2; stop=sr/2, length=length(ts)), # esto es aprox
fftshift(abs.(x2s_fft)))
Podríamos armar una funcioncita que grafica siguiendo estos pasos
function plotspec(x, sr, tmax; kw...)
ts = range(0; stop = tmax, step=1/sr)
xs = x.(ts)
xs_fft = fft(xs)
plot(
range(-sr/2; stop=sr/2, length=length(ts)), # esto es aprox
fftshift(abs.(xs_fft));
xlabel = "f (Hz)",
kw...)
end
plotspec(x2, sr, 3; title="|X2(f)|")