Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

how to extract frequency associated with fft values in python

Tags:

python

numpy

fft

I used fft function in numpy which resulted in a complex array. How to get the exact frequency values?

like image 905
ria Avatar asked Sep 12 '10 12:09

ria


People also ask

How do you find the frequency spectrum of a signal in Python?

The frequency components of the periodic signal are plotted in the horizontal axis and amplitude component of the periodic signal is plotted in the vertical axis. In Python, the magnitude_spectrum() method in the pyplot module of Python matplotlib library plots the magnitude spectrum of a periodic signal.

What is Fftfreq Python?

fftfreq() method, we can compute the fast fourier transformation frequency and return the transformed array by using this method. Syntax : scipy.fftfreq(n, freq) Return : Return the transformed array.


2 Answers

np.fft.fftfreq tells you the frequencies associated with the coefficients:

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) 

The OP asks how to find the frequency in Hertz. I believe the formula is frequency (Hz) = abs(fft_freq * frame_rate).

Here is some code that demonstrates that.

First, we make a wave file at 440 Hz:

import math import wave import struct  if __name__ == '__main__':     # http://stackoverflow.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() 

This creates the file test.wav. Now we read in the data, FFT it, find the coefficient with maximum power, and find the corresponding fft frequency, and then convert to 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 
like image 200
unutbu Avatar answered Sep 24 '22 23:09

unutbu


Frequencies associated with DFT values (in python)

By fft, Fast Fourier Transform, we understand a member of a large family of algorithms that enable the fast computation of the DFT, Discrete Fourier Transform, of an equisampled signal.

A DFT converts a list of N complex numbers to a list of N complex numbers, with the understanding that both lists are periodic with period N.

Here we deal with the numpy implementation of the fft.

In many cases, you think of

  • a signal x defined in the time domain of length N, sampled at a constant interval dt,
  • its DFT X (here specifically X = np.fft.fft(x)), whose elements are sampled on the frequency axis with a sample rate dw.

Some definition

  • the period (aka duration) of the signal x, sampled at dt with N samples is is

    T = dt*N 
  • the fundamental frequencies (in Hz and in rad/s) of X, your DFT are

    df = 1/T dw = 2*pi/T # =df*2*pi 
  • the top frequency is the Nyquist frequency

    ny = dw*N/2 

    (and it's not dw*N)

The frequencies associated with a particular element in the DFT

The frequencies corresponding to the elements in X = np.fft.fft(x) for a given index 0<=n<N can be computed as follows:

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

or in a single sweep

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

if you prefer to consider frequencies in Hz, s/w/f/

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

Using those frequencies

If you want to modify the original signal x -> y applying an operator in the frequency domain in the form of a function of frequency only, the way to go is computing the w's and

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

Introducing np.fft.fftfreq

Of course numpy has a convenience function np.fft.fftfreq that returns dimensionless frequencies rather than dimensional ones but it's as easy as

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

Because df = 1/T and T = N/sps (sps being the number of samples per second) one can also write

f = np.fft.fftfreq(N)*sps 
like image 34
gboffi Avatar answered Sep 24 '22 23:09

gboffi