Signal processing with Fourier transform

Question

I need to process a periodic signal and obtain its frequency. This can be easily done with scipy.fft and it works fine.

However, I have a particular signal that is not strictly periodic. The period changes slightly but randomly. Applying the fft to this signal is quite hard. Many peaks are obtained and I cannot extrapolate only the range of frequencies that are near the (1/period) I am interested in.

How can I do this?

This is a simple code snippet of what I am doing:

df = pd.read_csv('data.txt', header=None)
x = np.asarray(df.iloc[:, 0])
y = np.asarray(df.iloc[:, 1])

yf = fft(y)
xf = fftfreq(len(y))

plt.plot(xf, abs(yf))

You can find an example of such signal at the following GitHub repository, inside the README file: https://github.com/crazydecibel/Stack-Overflow-question

Answer 1

What about taking the weighted average of frequency of top energy bins?

import numpy as np
import matplotlib.pyplot as plt

your_data = '15.042,...' # your github data
t, y = np.array([list(map(float, row.split(','))) for row in your_data.split()]).T

# here is the solution in terms of t and y.
Y = np.fft.fftshift(np.fft.fft(y))
Ts = (t[-1] - t[0]) / (len(t) + 1) # sampling period
Fs = 1 / Ts # sampling frequency
f = np.fft.fftshift(np.fft.fftfreq(len(y))) * Fs
# get top 5%
top = np.argsort(abs(Y))[-10:]

plt.subplot(211)
plt.stem(f, abs(Y), 'o')
plt.stem(f[top], abs(Y[top]), '+r')
plt.xlim(-0.05 * Fs, 0.05 * Fs)

fc = np.sum(abs(f[top]*Y[top]**2))/np.sum(abs(Y[top])**2)
plt.subplot(212)
plt.plot(t, y)
plt.plot(t, np.sin(t*(2*np.pi*fc)))