Cedgar
Cedgar

Reputation: 33

Discrete Fourier Transform not working/very inefficient in python

I am trying to write a Discrete Fourier Transform function in python that will give the energy spectral density of a signal as an array (which I will then output graphically). I am doing this using matrix multiplications. My code seems to work for a small set of the data but takes a long time to process; for a large set of data e.g. a wav file it never completes the task. The function is currently:

from scipy.io import wavfile
import numpy as np
import matplotlib.pyplot as plt

def esd(data, fs):

a=[]


for i in range(int(np.size(data)/100)):

    dt = 1/(fs)

    fvec = np.arange(100*i , 100+(100*i) , 1)  
    nvec = np.arange(0 , np.size(data) , 1)

    Wconst = np.exp(-1j*2*np.pi/np.size(data))

    ematrix = Wconst**(np.outer(fvec,nvec))

    b = (dt**2)*(abs(np.inner(ematrix , data))**2)

    a = np.concatenate((a,b))




return a

Where data is the data set and fs is the sampling frequency. Is there a reason why this function is so slow or inefficient? It is designed to process the signal in 100 blocks of frequencies as otherwise the matrix would be extremely large.

Upvotes: 3

Views: 1018

Answers (1)

bnaecker
bnaecker

Reputation: 6430

This algorithm implements the "naive" discrete Fourier transform (DFT), by computing the Vandermonde frequency matrix. The problem is right here:

b = (dt ** 2) * abs(np.inner(ematrix, data)) ** 2

This naive implementation uses a direct matrix-vector multiplication, and has a running time of O(N**2), where N == data.size. So it'll do much worse as you get larger data, and probably not complete in a reasonable time for your WAV-file.

This is the whole point of using the Fast Fourier Transform algorithm, which takes advantage of lots of structure in the problem. Most notably, the FFT breaks the problem down recursively into smaller problems of size N / 2. That recursive structure is what gives the FFT a running time of O(N log N).

Upvotes: 3

Related Questions