在 python 中使用 FFT 重建原始信号

Question

我有一个分析生成的光谱，其中 x 轴代表 angular 频率，y 代表强度。频谱以某个频率值为中心，该频率值通常称为信号的中心频率（图中的蓝色图表）。我想对数据执行 IFFT 到时域，用高斯曲线切割它有用的部分，然后 FFT 回到原始域。我的问题是，在 IFFT(FFT(signal)) 之后，中心频率丢失了：我按形状取回了频谱，但它始终以 0 为中心（橙色图）。目前我对此的解决方案非常糟糕：我缓存原始 x 轴并在 FFT 调用时恢复它。这显然有很多缺点，我想改进它。下面我包含了一个演示问题的小演示。我的问题是：这可以用更优雅的方式解决吗？有没有办法在这个过程中不丢失中心频率？

import numpy as np
from scipy.fftpack import fft, ifft, fftshift, fftfreq
from scipy.interpolate import interp1d
import matplotlib.pyplot as plt

C_LIGHT = 299.793

def generate_data(start, stop, center, delay, GD=0, resolution=0.1):
    window = 8 * np.log(2) / 50
    lamend = 2 * np.pi * C_LIGHT / start
    lamstart = 2 * np.pi * C_LIGHT/stop
    lam = np.arange(lamstart, lamend + resolution, resolution) 
    omega = 2 * np.pi * C_LIGHT / lam 
    relom = omega - center
    _i = np.exp(-(relom) ** 2 / window)
    i = 2 * _i + 2 * np.cos(relom * GD + (omega * delay)) * np.sqrt(_i * _i)
    return omega, i


if __name__ == '__main__':

    # Generate data
    x, y = generate_data(1, 3, 2, 800, GD=0)

    # Linearly interpolate to be evenly spaced
    xs = np.linspace(x[0], x[-1], len(x))
    intp = interp1d(x, y, kind='linear')
    ys = intp(xs)
    x, y = xs, ys
    plt.plot(x, y, label='original')

    # IFFT 
    xt = fftfreq(len(x), d=(x[0]-x[1])/(2*np.pi))
    yt = ifft(y)
    # plt.plot(xt, np.abs(yt))

    # FFT back
    xf = fftshift(fftfreq(len(xt), d=(xt[0]-xt[1])/(2*np.pi)))
    yf = fft(yt)
    plt.plot(xf, np.abs(yf), label='after transforms')
    plt.legend()
    plt.grid()
    plt.show()

Answer 1

我认为 fftfreq 并不像您认为的那样。 fft(ifft(y) 的 xf 与 x 相同，您不应尝试重新计算它。转到另一个域然后再返回时，x 轴不会改变。

另外，请注意 fftfreq returns 给定长度和给定样本间距的信号的离散傅立叶变换的频域坐标。它不会做相反的事情，在应用逆离散傅立叶变换后，您不能用它来确定空间域中的坐标。（它的间距returns是有效的，但坐标组不是。）

    plt.plot(x, y, label='original')

    # IFFT 
    yt = ifft(y)
    # plt.plot(np.abs(yt))

    # FFT back
    yf = fft(yt)
    plt.plot(x, np.real(yf), label='after transforms')
    plt.legend()
    plt.grid()
    plt.show()

您的代码的另一个问题是 ifft(y) 假定沿 x 轴的一组固定值。您的 x 与此不符。所以，你得到的空域信号是没有意义的。

运行你的代码，我看到 x 运行s 从 3.0 到 1.0，步长为 0.0004777。您必须扩充数据，使值运行从 0.0 到 6.0，区域 (3.0, 6.0) 是区域 (0.0, 3.0) 的共轭对称副本。根据频域的周期性（F[n]==F[n+N]，其中 N 是样本数），该区域对应于负频率。用零填充区域 (0.0, 1.0)。

给定频域中的标准化 x 轴，xf = fftfreq(len(xt), d=(xt[1]-xt[0])) 应该重建 x 轴。但是您需要适当地计算 xt：xt = np.linspace(0, 1/(x[1]-x[0]), len(x), endpoint=False)（使用 x 标准化 DFT 频率轴）。

在 python 中使用 FFT 重建原始信号

Reconstruct original signal with FFT in python

python

numpy

fft

scipy