I sometimes deal with time series and waveform data, so I tried some smoothing. I will write down the three methods together as a memorandum.
Generate the waveform data to be used this time
import numpy as np
import matplotlib.pyplot as plt
N = 1024 #The number of samples
dt = 0.001 #Sampling cycle[s]
f1, f2 = 30, 90 #frequency[Hz]
a1, a2 = 1.5, 0.8 #amplitude
t = np.arange(0, N*dt, dt) #time[s]
wave = a1*np.sin(2*np.pi*f1*t) + a2*np.sin(2*np.pi*f2*t) + 0.3 * np.random.randn(N) #signal
fig, ax = plt.subplots(figsize=(14.0, 6.0))
ax.plot(t, wave)
ax.set_xlabel("Time [s]")
ax.set_ylabel("Signal")
ax.grid()
plt.title("wave")
plt.show()
fig, ax = plt.subplots(figsize=(14.0, 6.0))
ax.plot(t, wave)
ax.set_xlim(0, 0.1)
ax.set_xlabel("Time [s]")
ax.set_ylabel("Signal")
ax.grid()
plt.title("enlarged view")
plt.show()
window = 5 #Moving average range
w = np.ones(window)/window
x = np.convolve(wave, w, mode='same')
fig, ax = plt.subplots(figsize=(14.0, 6.0))
ax.plot(t, x, label='moving average')
ax.plot(t, wave, alpha=0.3, label='wave')
ax.set_xlabel("Time [s]")
ax.set_ylabel("Signal")
ax.grid()
plt.legend()
plt.title("moving average window=5")
plt.show()
fig, ax = plt.subplots(figsize=(14.0, 6.0))
ax.plot(t, x, label='moving average')
ax.plot(t, wave, alpha=0.3, label='wave')
ax.set_xlim(0, 0.1)
ax.set_xlabel("Time [s]")
ax.set_ylabel("Signal")
ax.grid()
plt.legend()
plt.title("enlarged view")
plt.show()
import scipy.signal
#Degree of polynomial: 2
#Frame length: 5
x = scipy.signal.savgol_filter(wave, 5, 2, deriv=0)
fig, ax = plt.subplots(figsize=(14.0, 6.0))
ax.plot(t, x, label='Savitzky-Golay')
ax.plot(t, wave, alpha=0.3, label='wave')
ax.set_xlabel("Time [s]")
ax.set_ylabel("Signal")
ax.grid()
plt.legend()
plt.title("Savitzky-Golay deriv=0")
plt.show()
fig, ax = plt.subplots(figsize=(14.0, 6.0))
ax.plot(t, x, label='Savitzky-Golay')
ax.plot(t, wave, alpha=0.3, label='wave')
ax.set_xlim(0, 0.1)
ax.set_xlabel("Time [s]")
ax.set_ylabel("Signal")
ax.grid()
plt.legend()
plt.title("enlarged view")
plt.show()
x = scipy.signal.savgol_filter(wave, 5, 2, deriv=1)
fig, ax = plt.subplots(figsize=(14.0, 6.0))
ax.plot(t, x)
ax.set_xlabel("Time [s]")
ax.set_ylabel("Signal")
ax.grid()
plt.title("Savitzky-Golay deriv=1")
plt.show()
fig, ax = plt.subplots(figsize=(14.0, 6.0))
ax.plot(t, x)
ax.set_xlim(0, 0.1)
ax.set_xlabel("Time [s]")
ax.set_ylabel("Signal")
ax.grid()
plt.title("enlarged view")
plt.show()
x = scipy.signal.savgol_filter(wave, 5, 2, deriv=2)
fig, ax = plt.subplots(figsize=(14.0, 6.0))
ax.plot(t, x)
ax.set_xlabel("Time [s]")
ax.set_ylabel("Signal")
ax.grid()
plt.title("Savitzky-Golay deriv=2")
plt.show()
fig, ax = plt.subplots(figsize=(14.0, 6.0))
ax.plot(t, x)
ax.set_xlim(0, 0.1)
ax.set_xlabel("Time [s]")
ax.set_ylabel("Signal")
ax.grid()
plt.title("enlarged view")
plt.show()
3.1 Fast Fourier Transform
x = np.fft.fft(wave)
x = np.abs(x) #Convert complex numbers to absolute values
x = x / N * 2 #Amplitude adjustment
fq = np.linspace(0, 1.0/dt, N) #Frequency adjustment
fig, ax = plt.subplots(figsize=(14.0, 6.0))
ax.plot(fq[: int(N / 2)], x[: int(N / 2)])
ax.set_xlabel("Frequency [Hz]")
ax.set_ylabel("Swing")
ax.grid()
plt.title("Fast Fourier Transform")
plt.show()
#Frequency 30,90, amplitude 1.5, 0.Peak is seen around 8
3.2 Inverse Fast Fourier Transform without Noise
threshold = 0.6 #Amplitude threshold
x = np.fft.fft(wave)
x_abs = np.abs(x)
x_abs = x_abs / N * 2
x[x_abs < threshold] = 0
x = np.fft.ifft(x)
x = x.real #Extract only the real part from a complex number
fig, ax = plt.subplots(figsize=(14.0, 6.0))
ax.plot(t, x, label='IFFT')
ax.plot(t, wave, alpha=0.3, label='wave')
ax.set_xlabel("Time [s]")
ax.set_ylabel("Signal")
ax.grid()
plt.legend()
plt.title("Inverse Fast Fourier Transform")
plt.show()
fig, ax = plt.subplots(figsize=(14.0, 6.0))
ax.plot(t, x, label='IFFT')
ax.plot(t, wave, alpha=0.3, label='wave')
ax.set_xlim(0, 0.1)
ax.set_xlabel("Time [s]")
ax.set_ylabel("Signal")
ax.grid()
plt.legend()
plt.title("enlarged view")
plt.show()
When I have time, I would like to add the theory part and describe other methods.