import pandas as pd
import numpy as np
import matplotlib.dates as mdates
from datetime import datetime
import matplotlib.pyplot as plt
df = pd.read_csv('datasets/ocean-temperature-westfjords-2024/cmems_obs-ins_nws_phybgcwav_mynrt_na_irr_1763979809373_noheader.csv')
time = df["time"]
temperature = df["value"]
iso_date_string = "2023-05-29T10:30:00Z"
parsed_date = datetime.fromisoformat(iso_date_string[:-1]) # Removing the 'Z' at the end
parsed_time = [0] * len(time)
for i in range(len(time)):
parsed_time[i] = datetime.fromisoformat(time[i])
I continued to use Neil's FFT code, but I'm not sure what to make the rate.
rate = temperature.size
fft = np.fft.fft(temperature)
frequencies = np.fft.fftfreq(len(temperature),d=1/rate)
plt.plot(frequencies,abs(fft))
plt.xlim(0,1)
plt.title('FFT')
plt.xlabel('frequency (Hz)')
plt.ylabel('magnitude')
plt.show()
I'm not sure what I need to do to make this work.
Wave height dataset¶
However, I found a nice example of ocean wave height data that is probably more amenable to the Fourier transform. It was measured every 0.5 seconds with a buoy on a Norwegian research vessel. This vessel is very similar to the research vessels my father worked on his whole career as a fisheries expert at the Icelandic Marine & Freshwater Research Institute.
import numpy as np
import matplotlib.pyplot as plt
#read and plot wave elevation data
elevation_dataset = np.loadtxt('datasets/wave-height/elevation_vs_time.txt')
time = elevation_dataset[:,0]
elevation = elevation_dataset[:,1]
n_points = time.size
plt.plot(time,elevation)
plt.show()
[<matplotlib.lines.Line2D at 0xec2611017890>]
Let's plot a shorter section of time to see the waveform:
plt.plot(time,elevation)
plt.xlim(0,200)
plt.title('Wave elevation')
plt.xlabel('Time [s]')
plt.ylabel('Elevation [m]')
plt.show()
Now it's time to transform the time domain data into the frequency domain. This is Python code from Milan Stanko, based on the document "Frequency Domain Using Excel" by Larry Klingenberg at San Fransisco State University.
import scipy.signal as signal
fft_peaks = signal.find_peaks(fft_amplitude, prominence=1)
print(fft_peaks)
(array([ 147, 179, 212, 3884, 3917, 3949]), {'prominences': array([1.20333831, 1.55364599, 1.05997408, 1.05997408, 1.55364599,
1.20333831]), 'left_bases': array([ 49, 49, 199, 1924, 1924, 3933]), 'right_bases': array([ 163, 1924, 1924, 3897, 4047, 4047])})
Finding the peaks doesn't quite work, but here is the FFT plot:
plt.plot(fft_frequency,fft_amplitude)
# plt.scatter(0,fft_peaks[0][0],color = 'red')
plt.xlim(-0.3,0.3)
plt.show()
The highest peak occurs at about 0.8 Hz. That's the most common wave frequency on this ocean voyage.