千字解析线性卷积、循环卷积与周期卷积

我们在对音频信号进行处理时很少会直接在时域进行滤波运算，一般都会变换到频域进行计算。究其原因，我想很多人都听说过：FFT可以快速实现卷积运算，我们今天就介绍下卷积和FFT的爱恨情仇。当然，我们这里的卷积指的是信号处理里面的卷积，深度学习里面的卷积不在今天讨论范围之内。

三种卷积

卷积是一种在信号处理和图像处理中广泛应用的运算方法，它可以对信号或图像进行特征提取和信息融合。在卷积运算中，线性卷积和循环（圆周)卷积与周期卷积是常见的卷积方式，它们在数学原理和应用场景上有着不同的特点。线性卷积和循环卷积在某些条件下是等效的，循环卷积是周期卷积的一个特例，下面我们分别介绍这三种卷积。

线性卷积是一种常见的卷积方式，它是通过在输入信号上滑动一个滤波器（也称为卷积核或卷积矩阵）来实现的。在每个位置，输入信号和滤波器的对应元素相乘，然后将乘积相加得到卷积结果，其计算公式为：

循环卷积是一种在离散周期信号处理中常用的卷积方式，它考虑了信号在周期边界上的周期性。在圆周卷积中，输入信号和滤波器的数列在进行卷积时可以循环使用，从而避免了线性卷积中输入信号和滤波器在边界处的填充问题，其计算公式为：

注意到循环卷积n的范围不再是-∞到+∞，并且计算公式后面多了模运算R。

周期卷积是对具有周期性的信号执行的卷积运算，其计算公式与(2)相同，但是它要求两个参与计算的信号都是周期性的并且周期一样，如果不满足则对它们进行周期延拓，并且周期卷积只在一个周期内进行计算。

光看公式有点抽象，我们举个例子，假设离散序列x[n]，如下所示

我们还有一个单位延迟序列h[n]，如下所示

那么，x[n]与h[n]线性卷积结果如下，相当于对x[n]向右进行平移了一个单位 。

x[n]与h[n]循环卷积结果如下，整体上是向右平移了一个单位，但是本该在序列最右边的元素因为循环移位出现在了最左边的位置。

最后计算x[n]与h[n]周期卷积结果，因为周期卷积要求信号是周期的，我们先对x[n]和h[n]进行周期延拓，然后在一个周期内计算它们的循环卷积，如下所示。可以看出

• 周期卷积是线性卷积以一定序列长度为周期的周期延拓
• 周期卷积取主值序列即为循环卷积

卷积与混叠

假设卷积序列x[n]与h[n]的长度分别为M和N，卷积后的结果y[n]长度为L，那么线性卷积结果长度为N+M-1。在计算周期卷积时，需要对两个有限长序列的线性卷积的结果按照卷积周期L进行延拓，当周期长度小于线性卷积结果的长度时，周期卷积序列进行周期延拓时必然有M+N-1-L点产生混叠，因此不发生混叠的条件是L≥N+M-1，从而可以得出以下结论：

• 当L<N+M-1时，循环卷积是线性卷积长度为L的混叠
• 当L=N+M-1时，循环卷积=线性卷积
• 当L>N+M-1时，循环卷积是线性卷积末尾补L-(N+M-1)个零

卷积与FFT

有了以上的了解我们就可以使用FFT加速线性卷积运算了，整体过程的伪代码如下。时域进行线性卷积运算，其时间复杂度为O(N^2)。但是使用FFT加速运算，其时间复杂度降到了O(NlogN)。

Input: Input signal x(n) of length M
       Impulse response h(n) of length L
Output: Linear convolution result y(n) of length N (N >= M + L - 1)

1. Zero-pad x(n) and h(n) to length N:
   if M > N, zero-pad h(n) with (M + L - 1 - N) zeros
   else if L > N, zero-pad x(n) with (M + L - 1 - N) zeros

2. Compute DFT of x(n) and h(n):
   X(k) = FFT(x(n)), k = 0, 1, ..., N-1
   H(k) = FFT(h(n)), k = 0, 1, ..., N-1

3. Compute Y(k) = X(k) * H(k), element-wise multiplication of X(k) and H(k)

4. Compute IDFT of Y(k) to obtain y(n):
   y(n) = IFFT(Y(k)), n = 0, 1, ..., N-1

5. Take the first N samples of y(n) as the linear convolution result
   y(n) = y(n), n = 0, 1, ..., N-1

块卷积

在许多应用中，比如对语音波形进行滤波，输入数据的长度可能是无限的。在这种情况下，计算整个输入信号的DFT可能是不切实际的，在滤波之前获得所有输入样本也会导致较长的延迟。这里我们可以使用块卷积，即将输入信号分成长度为L的片段，然后使用DFT将每个片段与FIR进行卷积，并通过合成滤波后的片段来获得所需的线性卷积结果。下面介绍常用的两种方法。

重叠相加法（overlap-add method）将输入信号分段为长度为L的片段，并将每个片段与长度为P的FIR进行卷积。输入一个片段与FIR的线性卷积将得到长度为(L + P – 1)的序列y[n]。因此，我们可以使用长度为(L + P – 1)的DFT来计算卷积，从而确保不会出现时间混叠。在滤波后的片段中，非零点总共有(P – 1)个重叠的点，将这些重叠的点相加来获得最终的输出。

在重叠相加法中，在计算每个片段后，我们需要存储 y[n] 的 (P – 1) 个值，并等待下一个数据段来添加重叠点。在某些情况下，这可能不是理想的做法，我们可以使用另一种方法，即重叠保存法（overlap-save method）。我们看到在循环卷积中，并非所有点都受到时间混叠的影响。每个片段的前 (P – 1) 个点受到时间混叠的影响，但是我们仍有 L – (P – 1) = (L – P + 1) 个点等于线性卷积。因此，在每个输出序列中位于 0 ≤ n ≤ P – 2 的部分将被丢弃，而剩余的样本将被保存以构造最终的滤波输出。

参考文献：

[1]. https://www.zhihu.com/question/394657296

[2]. https://ww2.mathworks.cn/help/signal/ug/linear-and-circular-convolution.html

[3]. https://dezeming.top/wp-content/uploads/2022/04/%E7%A6%BB%E6%95%A3%E4%BF%A1%E5%8F%B7%E7%9A%84%E5%91%A8%E6%9C%9F%E5%8D%B7%E7%A7%AF%E4%B8%8E%E7%BA%BF%E6%80%A7%E5%8D%B7%E7%A7%AF.pdf

[4]. https://ocw.mit.edu/courses/6-341-discrete-time-signal-processing-fall-2005/6e5190ef6e0d66c78bfdce2be6ce7125_lec16.pdf

[5]. https://thewolfsound.com/circular-vs-linear-convolution-whats-the-difference/#periodic-convolution

[6]. Discrete-Time Signal Processing， Oppenheim, Schafer & Buck

作者：Ryuk | 来源：公众号——语音算法组