Кратковременное преобразование Фурье в питоне - PullRequest
Новая
       13

Кратковременное преобразование Фурье в питоне

5 голосов
/ 22 февраля 2012

Я хочу получить частоту с максимальной мощностью для каждого момента в файле WAV. Поэтому я написал STFT на Python, используя fft от scipy. Я использовал оконную функцию kaiser от scipy. Все выглядит отлично, но мой вывод выглядит странно. У него есть очень маленькие цифры, а некоторые очень высокие.

вот вывод для одного файла wav: http://pastebin.com/5Ryd2uXj и вот код в Python: 

import scipy, pylab
import wave
import struct
import sys

def stft(data, cp, do, hop):
    dos = int(do*cp)
    w = scipy.kaiser(dos,12) //12 is very high for kaiser window
    temp=[]
    wyn=[]
    for i in range(0, len(data)-dos, hop):
        temp=scipy.fft(w*data[i:i+dos])
        max=-1
        for j in range(0, len(temp),1):
            licz=temp[j].real**2+temp[j].imag**2
            if( licz>max ):
                max = licz
                maxj = j
        wyn.append(maxj)
    #wyn = scipy.array([scipy.fft(w*data[i:i+dos])
        #for i in range(0, len(data)-dos, 1)])
    return wyn

file = wave.open(sys.argv[1])
bity = file.readframes(file.getnframes())
data=struct.unpack('{n}h'.format(n=file.getnframes()), bity)
file.close()

cp=44100 #sampling frequency
do=0.05 #window size
hop = 5

wyn=stft(data,cp,do,hop)
print len(wyn)
for i in range(0, len(wyn), 1):
    print wyn[i]

1 Ответ

5 голосов
/ 22 февраля 2012

Фактический FT синусоидальной волны представляет собой пару дельта-функций, равноудаленных от 0-частоты.Для дискретной функции (выборки) это повторяется каждые fs (частота дискретизации) в частотной области.Небольшие ошибки в вычислении БПФ будут означать, что эти две дельты (FT вашей синусоидальной волны) не будут точно одинаковой высоты, поэтому ваш алгоритм просто выбирает более высокую.

Функция scipy FFT даст вам частотные компонентыс доменом [0, fs].Поскольку (как я упоминал выше) это периодическое действие, эти значения также можно переназначить как [-fs/2, fs/2], поменяв местами результат в центральной точке - для этого используйте fftshift.Похоже, вас могут интересовать только положительные частоты, поэтому вы можете просто отбросить вторую половину результата вашего БПФ.

Из примечаний scipy.fftpack.fft:

Упаковка результата является «стандартной»: если A = fft (a, n), то A [0] содержит член нулевой частоты, A [1: n / 2 + 1] содержит члены с положительной частотой, а A [n / 2 + 1:] содержит члены с отрицательной частотой в порядке убывания отрицательной частоты.Таким образом, для 8-точечного преобразования частоты результата равны [0, 1, 2, 3, 4, -3, -2, -1].

...