Manual fft not giving me same results as fft

2019-08-29 06:01发布

问题:

import numpy as np
import matplotlib.pyplot as pp

curve = np.genfromtxt('C:\Users\latel\Desktop\kool\Neuro\prax2\data\curve.csv',dtype =     'float', delimiter = ',')
curve_abs2 = np.empty_like(curve)
z = 1j
N = len(curve)
for i in range(0,N-1):
    curve_abs2[i] =0
    for k in range(0,N-1):
        curve_abs2[i] += (curve[i]*np.exp((-1)*z*(np.pi)*i*((k-1)/N)))

for i in range(0,N):
    curve_abs2[i] = abs(curve_abs2[i])/(2*len(curve_abs2))

#curve_abs = (np.abs(np.fft.fft(curve)))
#pp.plot(curve_abs)
pp.plot(curve_abs2)
pp.show()

The code behind # gives me 3 values. But this is just ... different

Wrong ^^ this code: http://www.upload.ee/image/3922681/Ex5problem.png

Correct using numpy.fft.fft(): http://www.upload.ee/image/3922682/Ex5numpyformulas.png

回答1:

There are several problems:

  1. You are assigning complex values to the elements of curve_abs2, so it should be declared to be complex, e.g. curve_abs2 = np.empty_like(curve, dtype=np.complex128). (And I would recommend using the name, say, curve_fft instead of curve_abs2.)

  2. In python, range(low, high) gives the sequence [low, low + 1, ..., high - 2, high - 1], so instead of range(0, N - 1), you must use range(0, N) (which can be simplified to range(N), if you want).

  3. You are missing a factor of 2 in your formula. You could fix this by using z = 2j.

  4. In the expression that is being summed in the inner loop, you are indexing curve as curve[i], but this should be curve[k].

  5. Also in that expression, you don't need to subtract 1 from k, because the k loop ranges from 0 to N - 1.

  6. Because k and N are integers and you are using Python 2.7, the division in the expression (k-1)/N will be integer division, and you'll get 0 for all k. To fix this and the previous problem, you can change that term to k / float(N).

If you fix those issues, when the first double loop finishes, the array curve_abs2 (now a complex array) should match the result of np.fft.fft(curve). It won't be exactly the same, but the differences should be very small.

You could eliminate that double loop altogether using numpy vectorized calculations, but that is a topic for another question.