使用pyfftw模块计算高斯函数傅里叶变换异常结果求助

问题分析 #

你的代码存在几个关键问题，导致结果不符合预期：

• FFT结果为复数：pyfftw.interfaces.numpy_fft.fft返回复数数组，直接绘图会显示无意义内容。
• 频率轴错误：直接用原空间域的x轴代替频率域的k轴，完全不符合傅里叶变换的频率刻度逻辑。
• 未做频谱移位：FFT默认输出是从0频率开始的，需要移位到中心才能和理论上的对称高斯函数对齐。
• 幅度未匹配理论值：离散FFT的输出幅度需要结合采样参数调整，才能和连续傅里叶变换的解析解匹配。
• 绘图函数调用错误：plt.show缺少括号，无法触发图像显示。

修正后的代码 #

import numpy as np
import pyfftw
import matplotlib.pyplot as plt
from math import sqrt, pi

c = 3
A = 3

# 定义采样参数
x_start = -50
x_end = 50
N = 100
x = np.linspace(x_start, x_end, N)
dx = x[1] - x[0]  # 采样间隔

# 定义原高斯函数
y = A * np.exp(-(x ** 2) / (2 * c ** 2))

# 生成正确的频率轴，计算理论傅里叶变换解析解
k = np.linspace(-np.pi/dx, np.pi/dx, N)
fr = A * sqrt(2 * pi) * c * np.exp(-(k ** 2) * c ** 2 / 2)

# 用pyfftw计算FFT并处理结果
fourier_transform = pyfftw.interfaces.numpy_fft.fft(y)
# 移位频谱到中心、取幅度、调整归一化系数
ft_shifted = np.abs(pyfftw.interfaces.numpy_fft.fftshift(fourier_transform)) * dx / (2 * pi)

# 绘图
plt.plot(x, y, color='b', label="原高斯函数")
plt.plot(k, ft_shifted, color='r', label="FFT计算结果")
plt.plot(k, fr, color='g', label="理论傅里叶变换")
plt.xlabel('x (空间域) / k (频率域)')
plt.ylabel('函数值')
plt.title("高斯函数的傅里叶变换")
plt.grid(True)
plt.legend()
plt.show()

关键修改说明 #

• 正确生成频率轴：根据采样间隔dx计算频率范围[-π/dx, π/dx]，确保频率刻度和理论解匹配。
• FFT结果处理：
  • fftshift：将FFT输出的频谱从0频率中心移到对称中心，和理论高斯的对称性对齐。
  • np.abs：取复数结果的模（原函数是实偶函数，FT理论上是实函数，也可使用np.real，取模更稳妥）。
  • 幅度归一化：乘以dx/(2π)，完成离散傅里叶变换到连续傅里叶变换的幅度匹配。
• 修复绘图触发：补上plt.show()的括号，确保图像正常显示。

内容的提问来源于stack exchange，提问作者Ghilas