使用 numpy fft 提取相位信息答案

【问题标题】：extracting phase information using numpy fft使用 numpy fft 提取相位信息
【发布时间】：2016-07-01 08:11:10
【问题描述】：

我正在尝试使用快速傅立叶变换来提取单个正弦函数的相移。我知道在纸面上，如果我们将函数的变换表示为 T，那么我们有以下关系：

但是，我发现虽然我能够准确地捕获余弦波的频率，但除非我以极高的速率进行采样，否则相位是不准确的。例如：

import numpy as np
import pylab as pl

num_t = 100000
t = np.linspace(0,1,num_t)
dt = 1.0/num_t
w = 2.0*np.pi*30.0
phase = np.pi/2.0

amp = np.fft.rfft(np.cos(w*t+phase))
freqs = np.fft.rfftfreq(t.shape[-1],dt)

print (np.arctan2(amp.imag,amp.real))[30]

pl.subplot(211)
pl.plot(freqs[:60],np.sqrt(amp.real**2+amp.imag**2)[:60])
pl.subplot(212)
pl.plot(freqs[:60],(np.arctan2(amp.imag,amp.real))[:60])
pl.show()

使用 num=100000 点，我得到的相位为 1.57173880459。

使用 num=10000 点，我得到的相位为 1.58022110476。

使用 num=1000 点，我得到的相位为 1.6650441064。

怎么了？即使有 1000 分，我每个周期也有 33 分，这应该足以解决它。有没有办法增加计算频率点的数量？有没有办法用“低”点数做到这一点？

编辑：从进一步的实验来看，我似乎每个周期需要约 1000 个点才能准确地提取相位。为什么？！

编辑 2：进一步的实验表明，准确性与每个周期的点数有关，而不是与绝对数有关。增加每个周期的采样点数会使相位更准确，但如果信号频率和采样点数都增加相同的倍数，则精度保持不变。

【问题讨论】：

如果增加点数会怎样？ 1,000,000？一千万？相位是否逐渐接近 pi/2？
我不是专家，但是随着采样率的降低，FFT 频率和相位结果中的噪声（不可避免地）会增加。较低的采样率 - 这表现为频率分量幅度和偏移的扩散在主成分阶段。随着其他分量的增加，您还会看到结果中的主频率幅度减小。
是的，随着点数的增加，相位会变得任意准确。我应该补充一点，似乎重要的不是绝对的点数，而是每周期的点数使它更准确。以 1kHz 采样 30Hz 信号 10 秒与以 1kHz 采样 1 秒的精度相同，但以 1kHz 采样 3Hz 波 1 秒与以 10kHz 采样 30Hz 波 1 秒的精度相同。
当你添加量化时会发生什么，例如样本中有 8 位还是 12 位？
我不确定我是否理解。量化...位是什么意思？

标签： python numpy fft

【解决方案1】：

您的分数在区间内分布不均，最后的分数加倍：0 与1 相同。显然，您获得的分数越多，这一点就越不重要，但仍然会出现一些错误。您可以完全避免它，linspace 对此有一个标志。它还有一个标志可以直接将dt 与数组一起返回给您。

做

t, dt = np.linspace(0, 1, num_t, endpoint=False, retstep=True)

而不是

t = np.linspace(0,1,num_t)
dt = 1.0/num_t

然后它工作:)

【讨论】：

啊，所以我的 dt （稍微）错了 1/N 倍？它似乎确实有效......这有点令人尴尬 - 显然需要在 linspace 上进行 RTFM。谢谢！
不，你的 dt 不是问题，这只是为了方便。问题是，你会两次考虑一个点，在此结束的时期和从那里开始的下一个时期
但是我怎样才能得到准确的真实世界数据呢？如果我采用实时时间序列，我不能保证最后一点与第一个不对应，可以吗？

【解决方案2】：

仅当输入信号在 FFT 长度内完全是整数周期时，未旋转 FFT 结果箱中的相位值才是正确的。您的测试信号不是，因此 FFT 测量的部分与测试正弦曲线端点之间的信号不连续性的相位差有关。较高的采样率会产生与正弦曲线略有不同的最后端点，因此可能会产生较小的不连续性。

如果您想减少此 FFT 相位测量误差，请创建您的测试信号，以便您的测试相位以测试矢量（不是第一个样本）的精确中心（样本 N/2）为参考，然后执行fftshift 操作（旋转 N/2），以便在生成的长度为 N 的 FFT 输入向量中的第一个点和最后一个点之间没有信号不连续性。

【讨论】：

你能举一个简短的例子吗？在我拥有实时数据且不一定能控制周期性的实际示例中，如何提取有意义的相位值？
对于任意信号，只需将您的 FFT 相位结果引用到 FFT 矢量的中心，而不是任何端点。您可以通过使用 fftshift 预处理步骤或通过后处理来执行此操作：根据 bin 编号是奇数还是偶数，您可能必须翻转相位的符号（这是另一个域中的 fftshift） .
在无法控制 FFT 孔径周期性的“现实世界”中，您通常会将数据窗口化，因此肯定要测量驼峰中心附近的相位（窗口的功能）。没有信号可以测量端点处的相位，因为常见的窗口函数通常在那里接近于零。
你能推荐一个讨论这个的参考吗？我真的不明白测量窗口中间附近的相位或将相位引用到 FFT 向量的中心是什么意思？
很公平，稍后我会在有时间适当地制定它时发布它。谢谢。

【解决方案3】：

这段代码可能会有所帮助：

def reconstruct_ifft(data):
"""
In this function, we take in a signal, find its fft, retain the dominant modes and reconstruct the signal from that

Parameters
----------
data : Signal to do the fft, ifft

Returns 
-------
reconstructed_signal : the reconstructed signal

"""
N = data.size
yf = rfft(data)

amp_yf = np.abs(yf) #amplitude

yf = yf*(amp_yf>(THRESHOLD*np.amax(amp_yf)))

reconstructed_signal = irfft(yf)

return reconstructed_signal

0.01 是您希望保留的 fft 幅度的阈值。使 THRESHOLD 更大（大于 1 没有任何意义），将给出更少的模式并导致更高的均方根误差，但确保更高的频率选择性。（请为python代码调整TABS）

【讨论】：