I'm trying to visualize a signal and its frequency spectrum in Julia.
I found the FFTW package that provides the FFT and DSP for the frequencies.
Here is what I'm trying, with a sinusoidal signal:
using Plots
using FFTW
using DSP
# Number of points
N = 2^14 - 1
# Sample rate
fs = 1 / (1.1 * N)
# Start time
t0 = 0
tmax = t0 + N * fs
# time coordinate
t = [t0:fs:tmax;]
# signal
signal = sin.(2π * 60 * t) # sin (2π f t)
# Fourier Transform of it
F = fft(signal)
freqs = fftfreq(length(t), fs)
freqs = fftshift(freqs)
# plots
time_domain = plot(t, signal, title = "Signal")
freq_domain = plot(freqs, abs.(F), title = "Spectrum")
plot(time_domain, freq_domain, layout = 2)
savefig("Wave.pdf")
I expected to see a nice plot with a peak in the 60 Hz, but all I got was a weird result:
I'm ignoring the negative frequencies for now.
How should I do that in Julia?
What you call fs
in your code is not your sampling rate but the inverse of it: the sampling period.
The function fftfreq
takes the sampling rate as its second argument. Since what you give as the second argument is the sampling period, the frequencies returned by the function are incorrectly scaled by (1/(Ts^2))
.
I renamed fs
to Ts
and changed the second argument to fftfreq
to the sampling rate 1.0/Ts
. I think you also need to shift the result of fft
.
# Number of points
N = 2^14 - 1
# Sample period
Ts = 1 / (1.1 * N)
# Start time
t0 = 0
tmax = t0 + N * Ts
# time coordinate
t = t0:Ts:tmax
# signal
signal = sin.(2π * 60 .* t) # sin (2π f t)
# Fourier Transform of it
F = fft(signal) |> fftshift
freqs = fftfreq(length(t), 1.0/Ts) |> fftshift
# plots
time_domain = plot(t, signal, title = "Signal")
freq_domain = plot(freqs, abs.(F), title = "Spectrum", xlim=(-1000, +1000))
plot(time_domain, freq_domain, layout = 2)
savefig("Wave.pdf")
If you love us? You can donate to us via Paypal or buy me a coffee so we can maintain and grow! Thank you!
Donate Us With