не получается корректно отобразить частотный спектр (Фурье)

Question

Необходимо получить частоты присутствующих в сигнале гармоник. Изначально дан файл с самим сигналом и частота дискретизации 500 МГц. Сигнал выглядит следующим образом

Для того, чтобы получить частоты написал следующий код:

data = pd.read_csv('signal_1.txt', header=None)
from scipy.fft import fft, fftfreq
fs = 500 # Переведем частоту в МГц

N = data.size

"""Прямое БПФ"""
X = fft(data)
freq = fftfreq(N, 1./fs) # Частоты полученных сигналов


f = freq[:N//2] 
A = abs(X[:N//2]) # Амплитуда сигнала

# Частотный спектр
plt.subplot(4, 2, 2)
plt.plot(f, A)
plt.title('Спектр')
plt.xlabel('Частота, МГц')
plt.ylabel('Амплитуда')
plt.grid()
plt.plot(f, A)

По итогу выводится следующая картина

Однако, как мне сказали, мы должны получить всего лишь 4 частоты. Я только начинаю работать с сигналами, но просмотрев уже достаточно много примеров все равно не могу понять, где я ошибся. Был бы рад, если бы мне указали на мою ошибку

Таблицу прикладываю

Answer 1

Датафрейм data выглядит для функции fft двумерным массивом. Поэтому fft запускает процедуру разложения Фурье каждой строки отдельно. Но! каждая строка матрицы data состоит из одного значения, поэтому преобразование Фурье ничего с ним не делает, и fft(data) фактически возвращает копию data. Чтобы убедиться в этом, вычислите значение выражения fft(data)-data - результат окажется столбцом из (комплексных) нулей.

Для преобразования Фурье вам достаточно превратить data в одномерный массив. Например, после загрузки данных из файла сделайте так: data = data.to_numpy().flatten(). После такого преобразования ваш код нарисовал вот какой спектр:

Вуаля!

Рекомендации

Вместо pandas.read_csv я рекомендую использовать numpy.loadtxt - работает гораздо бодрее и на выходе получается массив numpy.
Если вас интересует только спектр, то для его построения можно воспользоваться функцией scipy.signal.periodogram - она строит спектр плотности энергии, то есть модуль fft от квадрата амплитуды волны. Периодограмма не содержит информации о фазе.

Соответственно, мой первый вариант был таким:

wave = np.loadtxt('signal_1.txt')
freqs, spectrum = scipy.signal.periodogram(wave, fs=500)
plt.plot(freqs, spectrum)

Если же вам нужно именно разложение Фурье, то есть амплитуды и фазы, то достаточно функции numpy.fft.rfft, которая оптимизирована для разложения вещественной последовательности в ряд Фурье. В частности, эта функция не порождает отрицательные частоты, поэтому не требуется шаманство с N//2.

rfft_data = np.abs(np.fft.rfft(wave))
rfft_freqs = np.fft.rfftfreq(len(wave), d=1.0/500)
plt.plot(rfft_freqs, rfft_data)

БЛОГ НА HUSL

не получается корректно отобразить частотный спектр (Фурье)

Ответы (1 шт):