为什么NumPy的fft函数返回不正确的相位信息？

Question

我比较了Matlab和numpy中的相位和幅度谱。我认为Matlab的工作是正确的，但是numpy计算出了正确的幅度谱，但相位谱很奇怪。我该如何更改Python代码以通过numpy进行正确的fft计算？

Matlab代码:

fs = 1e4;
dt = 1 / fs;
t = 0:dt:0.5;
F = 1e3;

y = cos(2*pi*F*t);
S = fftshift(fft(y) / length(y));
f_scale = linspace(-1, 1, length(y)) * (fs / 2);

a = abs(S);
phi = (angle(S));

subplot(2, 1, 1)
plot(f_scale, a)
title('amplitude')

subplot(2, 1, 2)
plot(f_scale, phi)
title('phase')

Python:

import numpy as np
import matplotlib.pyplot as plt

fs = 1e4
dt = 1 / fs
t = np.arange(0, 0.5, dt)
F = 1e3

y = np.cos(2*np.pi*F*t)
S = np.fft.fftshift(np.fft.fft(y) / y.shape[0])
f_scale = np.linspace(-1, 1, y.shape[0]) * (fs / 2)

a = np.abs(S)
phi = np.angle(S)

plt.subplot(2, 1, 1, title="amplitude")
plt.plot(f_scale, a)

plt.subplot(2, 1, 2, title="phase")
plt.plot(f_scale, phi)

plt.show()

matlab输出

numpy输出

- Антон Пименов

Answer 1

在理解np.arange时存在一个问题。它会在到达所需值之前停止一个dt之前（您传递的区间在右侧是开放的）。如果定义

t = np.arange(0, 0.5+dt, dt)

一切都会正常工作。

- user8408080

Answer 2

正如另一篇答案所指出的那样，为了使Python绘图与matlab输出相匹配，您需要调整t数组，以使其具有与matlab代码中的t数组相同的值。

但是，如果您的目的是在信号中拥有一个整数周期，以便FFT仅具有两个非零值（在±输入频率处），则Python代码是正确的。 Python代码中的相位看起来很奇怪，因为除了与信号频率相关的所有傅里叶系数都是（理论上）0。 由于有限精度算术，系数最终变成数值“噪音”，振幅非常小且相位基本上是随机的。

- Warren Weckesser

---

感谢您详细的回答，现在我清楚地理解了为什么我有不同的返回结果。 - Антон Пименов

感谢您关于“小幅噪声”的回答，现在我通过以下方式过滤小幅振荡得到了完美的图形：fft_y = fft(y);abs_y=np.abs(fft_y);fft_y = np.where(abs_y<10,0,fft_y)。 - Joey Gao