Это правильный способ чтения FFT аудио файла?(питон + wav) - PullRequest
0 голосов
/ 04 февраля 2019

Аудиофайл представляет собой 16-битный монофонический аудиофайл PCM с различными частотами дискретизации и длиной 10–30 мс.

import struct
from pydub import AudioSegment
import numpy as np
import matplotlib.pyplot as plt
import scipy.fftpack

sound = AudioSegment.from_wav("3000hz.wav")

raw_data = sound.raw_data# needs to be mono
sample_rate = sound.frame_rate
sample_size = sound.sample_width
channels = sound.channels

fmt = "%ih" % sound.frame_count() * channels
amplitudes= struct.unpack(fmt, raw_data)
yVals = scipy.fftpack.fft(amplitudes)

plt.plot(abs(yVals[:(len(yVals)/2)-1]),'r')
plt.show()

Выход с файлом WAV с частотой 3000 Гц (взятым из онлайн-генератора синусоидальных волн) приводит кприлично выглядящий БПФ, но с 9000, а не 3000. Это в 3 раза больше, чем в других тестах.Это нормально?И правильный ли код?

1 Ответ

0 голосов
/ 06 февраля 2019

Вызывая plt.plot() только с массивом y и без соответствующего массива x, он будет использовать 0, 1, ..., N-1 в качестве значений x.Это не то, что мы на самом деле хотим, нам нужна частота на оси X.

Давайте обозначим значение x, которое вы видите на графике прямо сейчас, «индексом бина».Пусть длина массива будет N, а частота дискретизации будет fs.При расчете БПФ индекс бина 0 соответствует частоте 0 Гц.Следующий индекс бина 1 соответствует частоте fs / N Гц.Это связано с тем, что БПФ будет иметь значения N и перейдет от 0 Гц до fs Гц, поэтому каждый шаг равен fs / N Гц.Таким образом, следующая ячейка соответствует 2 * fs / N Гц и так далее.И последний интервал N-1 равен (N-1)/N * fs Гц, то есть почти fs Гц.

Если мы хотим создать график, в котором спектр амплитуды зависит от частоты, то нам нужно вручную создать частотувектор, который содержит реальную частоту для каждого индекса бина.К счастью, scipy.fftpack содержит функцию для этого: fftfreq:

freq = scipy.fftpack.fftfreq(n=N, d=1.0 / fs)

Затем мы можем изменить вызов на plt.plot(), чтобы использовать freq в качестве xзначения вместо 0 ... N-1:

plt.plot(freq, abs(yVals), 'r')

При этом пик должен быть в правильном положении.

Если вы хотите видеть только односторонний спектр, то вы можете обрезатьи freq и yVals, как вы уже сделали в коде в вопросе.

...