扩展IFFT

Question

假设我有一个波形 x，比如说：

x = [math.sin(W*t + Ph) for t in range(16)]

这个波形的 W 和 Ph 是随意的，然后我用下面的方法计算它的（实数）快速傅里叶变换（FFT） f：

f = numpy.fft.rfft(x)

我可以用以下方式得到原始的 x：

numpy.fft.irfft(f)

现在，如果我想把恢复出来的波形的范围向左和向右扩展一些样本，该怎么做呢？也就是说，我想得到一个波形 y，使得 len(y) == 48，并且 y[16:32] == x，而 y[0:16] 和 y[32:48] 是原始波形的周期性延伸。

换句话说，如果FFT假设它的输入是一个无限的函数 f(t)，在 t = 0, 1, ... N-1 的范围内取样，那么我该如何恢复 f(t) 在 t<0 和 t>=N 的值呢？

注意：我用一个完美的正弦波作为例子，但实际上 x 可以是任何东西：比如说任意信号，如 x = range(16) 或 x = np.random.rand(16)，或者是从随机的 .wav 文件中截取的任意长度的片段。

Answer 1

对于在FFT窗口或长度上是周期性的静态波形，你可以简单地把这个波形重复几遍，或者把经过IFFT(FFT())处理后的波形重复，来在时间上延长它。对于那些时间上被截断的波形，如果它们的来源在FFT窗口或长度上不是周期性的，FFT的结果就会是一个频谱和一个Sinc函数的卷积。所以，要恢复原始的未截断的频谱内容，就需要某种类似于去卷积的处理。由于这个去卷积过程很难或者根本不可能实现，所以通常会使用一种分析/再合成的方法，比如相位声码器处理或其他频率估计方法。然后，这些估计出来的频率可能和单个原始FFT结果中的频率不一样，可以把它们输入到一个正弦波合成器组，或者通过相位修改的IFFT混合，或者其他再合成的方法，来创建一个更长的波形，同时保持大致相同的频谱内容。

Answer 2

以下示例可以帮助你更好地理解该怎么做：

>>> x1 = np.random.rand(4)
>>> x2 = np.concatenate((x1, x1))
>>> x3 = np.concatenate((x1, x1, x1))
>>> np.fft.rfft(x1)
array([ 2.30410617+0.j        , -0.89574460-0.26838271j, -0.26468792+0.j        ])
>>> np.fft.rfft(x2)
array([ 4.60821233+0.j        ,  0.00000000+0.j        ,
       -1.79148921-0.53676542j,  0.00000000+0.j        , -0.52937585+0.j        ])
>>> np.fft.rfft(x3)
array([ 6.91231850+0.j        ,  0.00000000+0.j        ,
        0.00000000+0.j        , -2.68723381-0.80514813j,
        0.00000000+0.j        ,  0.00000000+0.j        , -0.79406377+0.j        ])

当然，得到三个句点最简单的方法就是把三个逆傅里叶变换的结果在时间域里拼接起来：

np.concatenate((np.fft.irfft(f),) * 3)

但是如果你想或者必须在频率域里这样做，你可以尝试以下方法：

>>> a = np.arange(4)
>>> f = np.fft.rfft(a)
>>> n = 3
>>> ext_f = np.zeros(((len(f) - 1) * n + 1,), dtype=f.dtype)
>>> ext_f[::n] = f * n
>>> np.fft.irfft(ext_f)
array([ 0.,  1.,  2.,  3.,  0.,  1.,  2.,  3.,  0.,  1.,  2.,  3.])

Answer 3

假设我需要把恢复的波形向左和向右扩展一些样本，比如说有一个波形 y，它的长度是 48，y[16:32] 等于 x，而 y[0:16] 和 y[32:48] 是原始波形的周期性扩展。

这里的 周期性 扩展其实就是 x，因为它是 周期性 的扩展。

换句话说，如果快速傅里叶变换（FFT）假设它的输入是一个在 t = 0, 1, ... N-1 之间取样的无限函数 f(t)，那么我怎么才能恢复 t<0 和 t>=N 的 f(t) 值呢？

这个“N点FFT假设”意味着你的信号是周期性的，周期为 N。这是因为你信号的基频函数都是周期性的，前面的 N 个样本和后面的 N 个样本其实就是主 N 个样本的复制。

如果你允许 W 取任何值，那么你的输入正弦波就不会是周期为 N 的。但这并不妨碍 FFT 函数把它分解成多个周期性的正弦波的和。而这些周期性正弦波的和也会有周期性 N。

显然，你需要重新考虑这个问题。

也许你可以利用线性预测。根据你片段的窗口自相关计算几个线性预测系数，然后用这些预测系数进行外推。不过，对于一个稳定的预测滤波器，预测会收敛到零，收敛的速度取决于你信号的类型。例如，对于白噪声，完美的线性预测系数都是零。在这种情况下，你只能向左和向右“外推”零。但对此你也无能为力。如果你有白噪声，片段中根本没有关于周围样本的信息，因为所有样本都是独立的（这就是白噪声的特点）。

这种线性预测实际上能够完美预测正弦波样本。所以，如果你的输入是 sin(W*t+p)，对于任意的 W 和 p，你只需要二阶的线性预测。对于更复杂的信号，我建议使用十阶或十六阶的线性预测。