如何确定用 scipy.fft 处理的信号分量波的实际相位角？

Question

我使用

scipy.fft

来转换信号并找到其分量。信号是随机选择的，所使用的代码如下所示。当手动构建分量波时，我对需要添加到每个分量的相位角感到困惑。

tt = 10*np.pi
N = 1000
fs = N/tt
space = np.linspace(0,tt,N)
signal = lambda t: 2*np.sin(t) + np.cos(2*t)

信号如下所示：

fft1 = fft.fft(signal(space))  
fft1t = fft1/len(space) # scaling the amplitudes
fft1t[1:] = fft1t[1:]*2 # scaling the amplitudes
freq = fft.fftfreq(space.shape[-1], d=1/fs)

茎图显示了正确识别的项及其幅度，并且从

freq

中识别出相应的频率，第 5 项为 0.1592，第 10 项为 0.3183。

plt.stem(np.abs(fft1t)[:int(len(fft1t)/2)]) # plotting only up to the Nyquist freq.
plt.xlim(-.5,20)
plt.grid()

现在，当我尝试手动重建信号并且绘图位于原始信号之上时，我得到的相位失配是根据 fft 结果计算出的相位角的 2 倍。

uk = space**0 * fft1t[0].real
for n in [5,10]: # 5th and 10th term of the Fourier series
    uk += (fft1t[n].real * np.cos(n * np.pi * space / (N/fs/2) + 0*np.angle(fft1t[n])) 
         + fft1t[n].imag * np.sin(n * np.pi * space / (N/fs/2) + 0*np.angle(fft1t[n])))
sns.lineplot(x=space, y=signal(space))
sns.lineplot(x=space, y=uk) 

当我添加一个相位角时，我得到：

uk = space**0 * fft1t[0].real
for n in [5,10]: # 5th and 10th term of the Fourier series
        uk += (fft1t[n].real * np.cos(n * np.pi * space / (N/fs/2) + 1*np.angle(fft1t[n])) 
             + fft1t[n].imag * np.sin(n * np.pi * space / (N/fs/2) + 1*np.angle(fft1t[n])))
sns.lineplot(x=space, y=signal(space))
sns.lineplot(x=space, y=uk)

当我添加两倍的相位角时，我得到了完美的匹配：

uk = space**0 * fft1t[0].real
for n in [5,10]: # 5th and 10th term of the Fourier series
        uk += (fft1t[n].real * np.cos(n * np.pi * space / (N/fs/2) + 2*np.angle(fft1t[n])) 
             + fft1t[n].imag * np.sin(n * np.pi * space / (N/fs/2) + 2*np.angle(fft1t[n])))
sns.lineplot(x=space, y=signal(space))
sns.lineplot(x=space, y=uk)

请帮助我理解为什么我需要将计算出的相位角加倍，以及确定要添加的相位角的一般规则是什么。在我最近发布的另一个问题（示例）中，我只需要添加一倍的相位角。这是否与采样空间相关，而不是从 [-T,T] 采样（如链接中的示例），我们在 [0,T] 中采样（如本示例中所示）？

Answer 1

如果您提供人们可以运行的完整代码，您将获得更多答案。

请注意，您使用的 numpy.linspace() 会在末尾生成一个您不想要的点。

您可以将信号重建为复指数的简单总和。您不应该对每个傅立叶系数应用不同的缩放比例。

我不明白为什么你需要任何相位角：

import numpy as np
from scipy import fft
from matplotlib import pyplot as plt

tt = 10*np.pi
N = 1000
fs = N/tt
space = np.linspace(0,tt,N,endpoint=False)       # Note that last point isn't included
space = space[0:N]
signal = lambda t: 2*np.sin(t) + np.cos(2*t)

# Compute FFT
fft1 = fft.fft(signal(space))  

# Invert FFT by direct calculation!
uk = 0
for n in range(0,N):
   uk += fft1[n] * np.exp((1.0j) * n * 2.0 * np.pi * space / tt )
uk = uk.real / len(space)

x = space
y = signal(space)
plt.plot(x,y)
plt.plot(x,uk,'x')
plt.show()