重叠-存储之卷积法 ( Overlap-save method, Overlap-discard method ) 是一种区块卷积 ( block convolution, sectioned convolution ),可以有效的计算一个很长的信号 xn 和一个 FIR 滤波器 hn 的离散卷积。
其中 hm 在 1, M 之外为零。
与重叠-相加之卷积法不同之处在于,重叠-存储之卷积法所算出的输出区块并不重叠 (因此计算上少了将输出区块相加所需的加法),而是每次用的输入区块有所重叠。因此实现时每次读取输入后需将和下一个输入重叠的部分存储起来,作为下一输入区块的开头部分,因此称为重叠-存储之卷积法。另外重叠-存储之卷积法也不需补零。
算法
概念上,这个做法是选用一个较短的适当长度 L 来切割 yn ,则因为 hn 是有限长度,因此在某一区块内的 yn 也只被有限长的 xn 区块(会比 yn 分区成的区块长一点)所决定。因此只要选择有影响的输入区块和 hn 卷积,再选择结果中适当的部分即可得到正确的输出区块。
则对于在 内的 n , 输出 yn 可写成
所以只需计算 n 在 中的 ykn + M - kL ,亦即 n 在 的 ykn 部分即可。因此每一段输出区块 ykn 的前 M-1 点可丢弃(discard)。
尽管一时看不出切割成区块的好处为何,但将 xkn 做 的周期延伸,
则 和 这两个卷积在 的部分相等。所以可以将线性卷积改以 点循环卷积计算,结果的 部分作为输出 yn 在 的部分。由于每段 xkn 原本就有 长,所以选择 的话输入 xn 就不需补零。
改以循环卷积计算后即可藉循环卷积定理
变换成三次 点快速傅里叶变换和 次乘法,使原本每段 O(N2) 的运算量减少至 O(N logN),速度大幅增加。
准代码
(Overlap-save algorithm for linear convolution)
//////// revised by fantastic ////////
N = length(x), M = length(h)
O = M – 1; // overlap length must be M-1
L = M; // >=1 is OK
P = O + L;
H = FFT(h, P); // just calc once
idx = - (O - 1); // starting index which is offset M-1 in matlab
while (idx <= N)
i1 = max(1, idx); // must be >= 1
i2 = min(N, idx+P-1); // must be <= N
yt = IFFT( FFT(x(i1:i2), P).*H, P );
y(idx:idx+P-M) = yt(M:P); // discard first M-1 values and concatenate the remaining
idx = idx + L;
end
y = y(1:M+N-1); // the first M+N-1 values are the convolution result
区块长度的选择
当 xn 的长度 N' 和 hn 的长度 M 相差太大时(例如 M < log2N' ),直接卷积(不透过循环卷积和 FFT )反而最快。而当 N' 和 M 差不多在同一个数量级时,不用分区,也就是只有一块长度 L = N' 的区块去做 FFT 即可。而当 N' 比 M 大了不少,却没大太多时,区块长度 L 就需要选择。除了与 N' 和 M 相关以外,也要考虑当两者相除有余数时,剩下一小段的输入可能会造成浪费。
相关条目
- 离散卷积
- 循环卷积
- 快速傅里叶变换
- 重叠-相加之卷积法
参考文献
- Rabiner, Lawrence R.; Gold, Bernard. Theory and application of digital signal processing. Englewood Cliffs, N.J.: Prentice-Hall. 1975: pp 65–67. ISBN 0-13-914101-4.
- Helms, H., Fast Fourier transform method of computing difference equations and simulating filters, IEEE Transactions on Audio and Electroacoustics, 1967, 15(2): 85–90
外部链接
- DSP class Fall 2005 Lecture08[永久失效链接] at The University of Texas at Arlington