我有一个时域表达式
f = -1j*H(t) * exp(-(1j*a+b)*t)
它可以用known properties(H
是Heaviside阶跃函数)进行解析变换。这种FT操作的结果是
其中w
是频率。
现在我使用this article中的技巧在Python中对f
进行数值傅里叶变换,并确认我确实得到了相同的分析结果F
:
import numpy as np
import matplotlib.pyplot as plt
t = np.linspace(-10,10,1e4) # time
w = np.linspace(-10,10,1e4) # frequency
b = 0.1
a = 1
H = lambda x: 1*(x>0) # heaviside function
# function in time
f = -1j*H(t)*np.exp(-(1j*a+b)*t)
# function in frequency (analytical work)
F = (w-a-1j*b)/((w-a)**2+b**2)
hann = np.hanning(len(t)) # hanning window
# function in frequency (numerical work)
F2 = 2/len(t)*np.fft.fft(hann*f)
plt.figure()
plt.plot(w,F.real,'b',label='analytical')
plt.plot(w,F2.real,'b--',label='fft')
plt.xlabel(r'$\omega$')
plt.ylabel(r'Re($F(\omega)$)')
plt.legend(loc='best')
plt.figure()
plt.plot(w,F.imag,'g',label='analytical')
plt.plot(w,F2.imag,'g--',label='fft')
plt.xlabel(r'$\omega$')
plt.ylabel(r'Im($F(\omega)$)')
plt.legend(loc='best')
plt.show()
但是Python的FFT函数似乎给了我完全错误的东西。当绘制F
和{
编辑:以下是绘图。。。
在这些图中并不明显,但是如果放大w=-10
和10
区域,会有小的振荡,可能是由于fft
算法造成的。
FFT算法计算DFT,它在第一个样本上有原点(包括空间域和频域)。你需要移动你的信号(在应用汉宁窗口之后),使t=0是最左边的样本,在计算FFT之后,你必须进行逆移位。在
MATLAB有},它们实现了这两个移位。NumPy必须有类似的功能。在
ifftshift
和{代码的另一个问题是计算DFT,并将其绘制在您计算的
w
给定的位置,但与计算DFT的实际频率无关。在这是您的代码,翻译成MATLAB,并修复了正确计算}*的代码。我希望这是有用的。需要注意的一点是你的}的错误。形状相似,但
F2
和{F
与F2
不匹配,我相信这不是由于F2
中的错误,而是由于你计算{F
的缩放和镜像方式不同。在脚注:
*我也改变了颜色,MATLAB的绿色太浅了。
相关问题 更多 >
编程相关推荐