CODEWITHSUNDEEP
pythonnumpysignal-processingfftdft
Shantanu Shinde
Shantanu Shinde

Reputation: 1012

Numpy fft function giving output different from the dft calculated using formula

I am trying to implement dft in python. I am using the standard formula: enter image description here

Here is my code:

k = np.array([np.arange(-50, 50)])
fs, xn = wavfile.read('voice_recording.wav')
nbits = 16
max_nbits = float(2**(nbits-1))
xn = xn / (max_nbits + 1.0)
xn = np.expand_dims(xn[:,0], axis=1)
N = len(xn)
n = np.array([np.arange(0, N)])
Xk = np.sum(xn*np.exp(((-1j*2*math.pi)/N)*np.matmul(n.T, k)), axis=0)

Here, xn is an audio signal read from a .wav file (voice_recording.wav). The code for FFT is:

Xk1 = np.fft.fftshift(np.fft.fft(xn, n=100, axis=0))

But both results are totally different, even though they should be same. DFT plot:

enter image description here

And FFT plot:

enter image description here

What am I doing wrong?

Upvotes: 0

Views: 548

Answers (1)

Cris Luengo
Cris Luengo

Reputation: 60504

Without downloading your data file, I presume that it has more than 100 samples. If so, then

np.fft.fft(xn, n=100, axis=0)

cuts off the first 100 samples and computes the FFT on those. That is, it does not compute the same thing as your code.

When I use xn = np.random.randn(100), and run your code, then both Xk and Xk1 are identical up to 1e-13 or so. This indicates that your code is correct.

To compute only a subset of k values using the FFT algorithm, first compute the full transform, then discard the values you don't want. For example:

Xk1 = np.fft.fft(xn, axis=0)
Xk1 = np.fft.fftshift(Xk1)
Xk1 = Xk1[(N//2 - 50):(N//2 + 50)]

Upvotes: 1

Related Questions