幅值和相位(角度)的fft.fft()输出与设置的值不对应

2022-03-28 00:00:00 python fft signal-processing

问题描述

我设置了一定幅度、频率和相位的正弦波,并尝试恢复幅度和相位:

import numpy as np
import matplotlib.pyplot as plt

N = 1000  # Sample points     
T = 1 / 800         # Spacing

t = np.linspace(0.0, N*T, N) # Time
frequency = np.fft.fftfreq(t.size, d=T) # Normalized Fourier frequencies in spectrum.

f0 = 25              # Frequency of the sampled wave
phi = np.pi/6       # Phase
A = 50              # Amplitude
s = A * np.sin(2 * np.pi * f0 * t - phi) # Signal

S = np.fft.fft(s)   # Unnormalized FFT

fig, [ax1,ax2] = plt.subplots(nrows=2, ncols=1, figsize=(10, 5))
ax1.plot(t,s,'.-', label='time signal')  
ax2.plot(freq[0:N//2], 2/N * np.abs(S[0:N//2]), '.', label='amplitude spectrum')

plt.show()

index, = np.where(np.isclose(frequency, f0, atol=1/(T*N))) # Getting the normalized frequency close to f0 in Hz)
magnitude = np.abs(S[index[0]]) # Magnitude
phase = np.angle(S[index[0]])   # Phase
print(magnitude)
print(phase)
phi
#21785.02149316858
#-1.2093259641890741
#0.5235987755982988

现在振幅应该是50,而不是21785,相位pi/6=0.524,而不是-1.2.

我是否误解了输出,或者上面链接中提到的帖子上的答案?


解决方案

  • 您需要使用以下两个更改之一将FFT归一化1/N(我使用的是第二个更改):
    S = np.fft.fft(s)-->;S = 1/N*np.fft.fft(s)
    magnitude = np.abs(S[index[0]])-->;magnitude = 1/N*np.abs(S[index[0]])
  • 不要使用index, = np.where(np.isclose(frequency, f0, atol=1/(T*N))),FFT不准确,最大幅度可能 不在f0,而使用np.argmax(np.abs(S)),它将提供 您的信号峰值将非常接近f0
  • np.angle搞砸了(我认为这是pi,pi/2 arctan偏移之一 事情)只需使用np.arctan(np.real(x)/np.imag(x))
  • 手动操作即可
  • 使用更多的点(我使N更高)并使T更小以获得更高的精度
  • 由于DFT(离散傅立叶变换)是双面的,并且在负频和正频都有峰值信号,所以正侧的峰值将只是实际幅度的一半。对于FFT,您需要将除f=0以外的每个频率乘以2才能进行计算。i在magnitude = np.abs(S[index])*2/N中乘以2
N = 10000
T = 1/5000
...
index = np.argmax(np.abs(S))
magnitude = np.abs(S[index])*2/N
freq_max = frequency[index]
phase = np.arctan(np.imag(S[index])/np.real(S[index])) 
print(f"magnitude: {magnitude}, freq_max: {freq_max}, phase: {phase}") print(phi)

输出:magnitude: 49.996693276663564, freq_max: 25.0, phase: 0.5079341239733628

相关文章