Python傅里叶逆变换相位不符：如何从Sinv中还原theta？

原始代码 #

import numpy as np
import matplotlib.pyplot as plt
from scipy.fft import fft, ifft, fftshift

def FourierTrans(S, t):
    fs = 1 / (t[5] - t[4])
    y = fft(S, len(S)*10)
    n = len(y)
    fshift = np.arange(-n/2,n/2)*(fs/n)
    yshift = fftshift(y)
    return fshift, yshift


def InvFourierTrans(Amp, freq):
    fs = 1/(freq[6]-freq[5])
    signal = ifft(Amp, len(Amp)*10)
    n = len(signal)
    time = np.arange(-n/2, n/2)*(fs/n)
    signal = fftshift(signal)
    return time, signal

plt.figure()
t = np.arange(-2.5, 2.5, 0.001)
sigma=np.diff(t)[0]*10
theta = 10000*t**2
y = np.exp(1j*theta) * np.exp(-t**2 / (2*sigma**2))
xnot, Snot= FourierTrans(y, t)
xinv, Sinv= InvFourierTrans(Snot, xnot)
plt.plot(xinv,np.angle(Sinv),linewidth=1,linestyle='-')
plt.plot(xinv, abs(Sinv)/max(abs(Sinv)),linewidth=1,linestyle='-')

plt.legend()
plt.grid()
plt.show()

问题描述 #

原本预期np.angle(Sinv)与np.angle(y)一致，但实际并非如此。已知abs(y)与abs(Sinv)除x轴偏移外完全匹配，如何从Sinv中还原theta？

原因分析 #

1. 补零操作与轴计算误差：变换时的补零（len(S)*10）以及手动计算频率/时间轴的方式，导致逆变换后的时间轴与原始t轴未对齐，且引入了线性相位偏移。
2. FFT的相位特性：fft和ifft的循环特性，加上fftshift的移位操作，会带来额外的周期性相位偏移。

解决方法 #

1. 修正变换函数的轴计算 #

用np.fft.fftfreq替代手动计算轴，确保频率和时间轴的准确性：

def FourierTrans(S, t):
    dt = t[1] - t[0]
    fs = 1 / dt
    n_pad = len(S)*10
    y = fft(S, n_pad)
    f = np.fft.fftfreq(n_pad, dt)
    fshift = fftshift(f)
    yshift = fftshift(y)
    return fshift, yshift

def InvFourierTrans(Amp, freq):
    df = freq[1] - freq[0]
    n_pad = len(Amp)
    signal = ifft(Amp, n_pad)
    t = np.fft.fftfreq(n_pad, df)
    tshift = fftshift(t)
    signal = fftshift(signal)
    return tshift, signal

2. 对齐时间轴并去除线性相位偏移 #

逆变换后的相位包含线性项，需要先对齐时间轴，再通过拟合去除线性偏移，还原原始二次相位：

# 原始信号
t = np.arange(-2.5, 2.5, 0.001)
sigma = np.diff(t)[0]*10
theta = 10000*t**2
y = np.exp(1j*theta) * np.exp(-t**2 / (2*sigma**2))

# 变换与逆变换
xnot, Snot = FourierTrans(y, t)
xinv, Sinv = InvFourierTrans(Snot, xnot)

# 对齐时间轴，截取与原始t匹配的区间
mask = (xinv >= t.min()) & (xinv <= t.max())
aligned_t = xinv[mask]
aligned_Sinv = Sinv[mask]

# 相位unwrap并拟合线性项
phase_unwrap = np.unwrap(np.angle(aligned_Sinv))
from scipy.optimize import curve_fit
def linear_func(x, a, b):
    return a*x + b

popt, _ = curve_fit(linear_func, aligned_t, phase_unwrap)
linear_phase = linear_func(aligned_t, *popt)
restored_theta = phase_unwrap - linear_phase

3. 验证结果 #

绘图对比原始theta与还原后的相位：

plt.figure(figsize=(10,6))
plt.plot(t, theta, label='原始theta', linewidth=1.5)
plt.plot(aligned_t, restored_theta, label='还原的theta', linewidth=1.5, linestyle='--')
plt.legend()
plt.grid(True)
plt.show()

修正后，还原的restored_theta会与原始theta基本重合，实现相位的准确还原。

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