Функция numpy fft, дающая результат, отличный от значения dft, рассчитанного по формуле

Question

Я пытаюсь реализовать dft в python. Я использую стандартную формулу:

Вот мой код:

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)

Здесь xn - это аудиосигнал, считанный с .wavфайл ( voice_recording.wav ). Код для БПФ:

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

Но оба результата совершенно разные, даже если они должны быть одинаковыми. График DFT:

И участок FFT:

Чтоя делаю не так?

Answer 1

Не загружая ваш файл данных, я предполагаю, что он содержит более 100 образцов. Если это так, то

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

отсекает первые 100 выборок и вычисляет БПФ для них. То есть он не вычисляет то же самое, что и ваш код.

Когда я использую xn = np.random.randn(100) и запускаю ваш код, тогда оба значения Xk и Xk1 идентичны до 1e-13 или около того,Это означает, что ваш код правильный.

Чтобы вычислить только подмножество k значений с использованием алгоритма FFT, сначала вычислите полное преобразование, а затем отбросьте значения, которые вам не нужны. Например:

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