transformada - señales de audio en python

cómo extraer la frecuencia asociada con los valores fft en python

Question

Usé la función fft en numpy, que dio como resultado una matriz compleja. ¿Cómo obtener los valores de frecuencia exactos?

Answer 1

Frecuencias asociadas con valores DFT (en python)

Por fft , Fast Fourier Transform, entendemos a un miembro de una gran familia de algoritmos que permiten el cálculo rápido de la DFT, Transformada de Fourier Discreta, de una señal equivipada.

Una DFT convierte una lista de N números complejos en una lista de N números complejos, en el entendimiento de que ambas listas son periódicas con el período N.

Aquí tratamos con la implementación numpy de la fft .

En muchos casos, piensas en

• una señal x definida en el dominio de tiempo de longitud N , muestreada a un intervalo constante dt ,
• su DFT X (aquí específicamente X = np.fft.fft(x) ), cuyos elementos se muestrean en el eje de frecuencias con una frecuencia de muestreo dw .

Alguna definición

• el período (también conocido como duración) de la señal x , muestreado en dt con N muestras es

  T = dt*N

• las frecuencias fundamentales (en Hz y en rad / s) de X , su DFT son

  df = 1/T
  dw = 2*pi/T # =df*2*pi

• la frecuencia superior es la frecuencia de Nyquist

  ny = dw*N/2

  (y no es dw*N )

Las frecuencias asociadas con un elemento particular en el DFT

Las frecuencias correspondientes a los elementos en X = np.fft.fft(x) para un índice dado 0<=n<N se pueden calcular de la siguiente manera:

def rad_on_s(n, N, dw):
    return dw*n if n<N/2 else dw*(n-N)

o en un solo barrido

w = np.array([dw*nif n<N/2 else dw*(n-N) for n in range(N)])

si prefiere considerar las frecuencias en Hz, s/w/f/

f = np.array([df*n if n<N/2 else df*(n-N) for n in range(N)])

Usando esas frecuencias

Si desea modificar la señal original x -> y aplicando un operador en el dominio de la frecuencia en forma de una función de frecuencia solamente, el camino a seguir es calcular las w y

Y = X*f(w)
y = ifft(Y)

Presentando np.fft.fftfreq

Por supuesto, numpy tiene una función de conveniencia np.fft.fftfreq que devuelve frecuencias adimensionales en lugar de frecuencias dimensionales, pero es tan fácil como

f = np.fft.fftfreq(N)*N*df
w = np.fft.fftfreq(N)*N*dw

Answer 2

La frecuencia es solo el índice de la matriz. En el índice n , la frecuencia es 2 πn / la longitud del conjunto (radianes por unidad). Considerar:

>>> numpy.fft.fft([1,2,1,0,1,2,1,0])
array([ 8.+0.j,  0.+0.j,  0.-4.j,  0.+0.j,  0.+0.j,  0.+0.j,  0.+4.j,
        0.+0.j])

el resultado tiene valores distintos de cero en los índices 0, 2 y 6. Hay 8 elementos. Esto significa

       2πit/8 × 0       2πit/8 × 2       2πit/8 × 6
    8 e           - 4i e           + 4i e
y ~ ———————————————————————————————————————————————
                          8

Answer 3

np.fft.fftfreq te dice las frecuencias asociadas con los coeficientes:

import numpy as np

x = np.array([1,2,1,0,1,2,1,0])
w = np.fft.fft(x)
freqs = np.fft.fftfreq(len(x))

for coef,freq in zip(w,freqs):
    if coef:
        print(''{c:>6} * exp(2 pi i t * {f})''.format(c=coef,f=freq))

# (8+0j) * exp(2 pi i t * 0.0)
#    -4j * exp(2 pi i t * 0.25)
#     4j * exp(2 pi i t * -0.25)

El OP pregunta cómo encontrar la frecuencia en Hertz. Creo que la fórmula es frequency (Hz) = abs(fft_freq * frame_rate) .

Aquí hay un código que demuestra eso.

Primero, hacemos un archivo de onda a 440 Hz:

import math
import wave
import struct

if __name__ == ''__main__'':
    # http://.com/questions/3637350/how-to-write-stereo-wav-files-in-python
    # http://www.sonicspot.com/guide/wavefiles.html
    freq = 440.0
    data_size = 40000
    fname = "test.wav"
    frate = 11025.0
    amp = 64000.0
    nchannels = 1
    sampwidth = 2
    framerate = int(frate)
    nframes = data_size
    comptype = "NONE"
    compname = "not compressed"
    data = [math.sin(2 * math.pi * freq * (x / frate))
            for x in range(data_size)]
    wav_file = wave.open(fname, ''w'')
    wav_file.setparams(
        (nchannels, sampwidth, framerate, nframes, comptype, compname))
    for v in data:
        wav_file.writeframes(struct.pack(''h'', int(v * amp / 2)))
    wav_file.close()

Esto crea el archivo test.wav . Ahora leemos los datos, lo convertimos en FFT, buscamos el coeficiente con la potencia máxima, buscamos la frecuencia fft correspondiente y luego convertimos a Hertz:

import wave
import struct
import numpy as np

if __name__ == ''__main__'':
    data_size = 40000
    fname = "test.wav"
    frate = 11025.0
    wav_file = wave.open(fname, ''r'')
    data = wav_file.readframes(data_size)
    wav_file.close()
    data = struct.unpack(''{n}h''.format(n=data_size), data)
    data = np.array(data)

    w = np.fft.fft(data)
    freqs = np.fft.fftfreq(len(w))
    print(freqs.min(), freqs.max())
    # (-0.5, 0.499975)

    # Find the peak in the coefficients
    idx = np.argmax(np.abs(w))
    freq = freqs[idx]
    freq_in_hertz = abs(freq * frate)
    print(freq_in_hertz)
    # 439.8975