區塊摺積

區塊折積,又稱區塊卷積分段折積分段卷積分塊卷積等,是一種計算線性離散卷积的方法,在信号处理以及線性非時變系統領域的應用層面上有相當高的價值。

區塊卷積主要常使用於計算一固定离散信号與另一非固定離散訊號間的線性卷積,且非固定訊號的長度通常較固定方長上許多,一個很具體的例子就是計算一訊號經有限長度濾波器後的輸出。
其主要具有兩個優點,第一,其計算複雜度較低,在輸入訊號長度為
N 時,區塊卷積的計算複雜度為\Omicron(N)\!,而直接對整段輸入訊號進行卷積計算的複雜度為\Omicron(N \log N)\!。
第二,當我們要以硬體進行實作時,區塊卷積僅需要單一固定的硬體架構,即可處理不同長度(甚至是無窮長)的輸入,可以像工廠的流水生產線一般連續的接受輸入並同時連續的吐出結果,
而如果要直接對整段輸入訊號進行卷積,則需要對不同的輸入長度設計不同的硬體運算架構,在應用上是不切實際並且沒有必要的。

目前對於區塊卷積,主流使用的英文專有詞彙為sectioned convolution,但在一些中文圈作者所撰寫的文章中,也可能會採用取「塊」字翻英後的block convolution。

計算
區塊卷積的計算核心概念,便是將較長的訊號進行分段,將原先的求卷積問題分解為規模較小的問題,再將結果進行整合,以得到原先問題的結果。

但根據將訊號分段以及將結果整合的方式,可以再將區塊卷積細分為兩個類別,可见下文说明。

卷積的性質觀察
首先,在對區塊卷積的計算法進行細分之前,我們需要先了解一些關於卷積的性質,以及關於切分後子問題的結果與母問題結果的關聯。

假設有一長度為 N 的離散訊號(通常很長):

x[n] ( n = 0 , 1 , \dots , (N-1) )

與一長度為 M 的離散訊號(通常很短):

h[n] ( n = 0 , 1 , \dots , (M-1) )

將兩者進行卷積後的結果 y[n] = x[n] * h[n],我們可以知道其長度應該為 (N+M-1) :

y[n] ( n = 0 , 1 , \dots , (N+M-2) )

根據卷積的定義,我們可以知道其中 y 的第 n' 點輸出 y[n'] ,
會受到相對應位置的 x[n'] ,以及其向前(M-1) 點的影響:

y[n'] = \sum_{m=0}^{M-1} x[n'- m] h[m] = x[n'- (M-1)]h[M-1]+ x[n' - (M-2)]h[M-2] + \dots + x[n']h[0]

這時如果我們將 x[n] 切出長度為 L 的小段,並假設起始位置為 n_0 :

x_L[n] = x[n] ( n = n_0 , n_0+1 , \dots , n_0+(L-1) )

我們可以發現,其與 h[n] 進行卷積後的結果 y_L[n] = x_L[n] * h[n],和原先母問題在對應區段的結果 y[n] 存在差異,
具體來說, y_L[n]的前段部分少考量了 x[n] 中 n部分的影響,例如在 y_L[n]第一點 y_L[n_0]
就僅考量一項 x 的影響,和 y[n_0] 相差許多:

y_L[n_0] = x[n_0]h[0] \neq y[n_0] = x[n_0- (M-1)]h[M-1]+ \dots + x[n_0]h[0]

同理,在 y_L 的後段部分也會有類似的問題,少考量了 x 在段落之後的影響,而只有計算到部分殘留影響的結果。
只有在 y_L 的中段部分,會與母問題結果 y 完全相同。
而要如何處理這個差異,並透過多段的 y_L 還原出 y ,便是以下兩種方法的最主要區別。

重疊-儲存之摺積法
重疊-儲存算法,便是捨棄子結果 y_L 的前後段,只取中段與母結果 y 相同的部分並進行拼接的演算法。

具體來說,重疊-儲存算法切分問題的主要考量為母結果 y ,
在實作這個算法的時候必須先決定將母結果 y 切分的長度。

假設決定將母結果切分為長度 L_y 的區間,
根據前面的討論,我們可以得知:

  1. x_L 的長度需要為 (L_y + (M-1)) 才能夠確保 L_y 長度的有效區段

2.在每次長度為 (L_y + 2(M-1)) 的計算結果 y_L 中,需要捨去前後各 (M-1) 長度的部分

3.在第一個區段,因為 x[n] 在 n的區域是沒有資訊的,所以 x_{L0} 的長度只需要 L_y 即可

綜合以上三點,我們可以寫出所需要的切分區段描述:

x_{L0}[n] = x[n] (n \in 0 \sim (L_y - 1))

x_{Lk}[n] = x[n] (n \in (kL_y - (M-1) ) \sim ((k+1)L_y - 1) , k \in 1 \sim (\left \lceil \frac{M+N-1}{L_y} \right \rceil - 1) )

以及結果的描述:

y[n] = y_{Lk}[n] (n \in kL_y \sim ((k+1)L_y - 1) , k \in 0 \sim (\left \lceil \frac{M+N-1}{L_y} \right \rceil - 1))

而重疊-儲存算法便是得名於其在輸入區段的切分上會有重疊的部分,須將前一個子問題的輸入後半先行儲存再作為下一個子問題的輸入前半。

這個算法的好處在於不需要額外的加法,並且如果在子問題求解時是採取直接計算卷積的方式,可以不用在訊號後面補零。

但缺點為同樣的分段數量下,每個子問題會需要處理較長的訊號,子問題計算量可能會較大一些。

重疊-相加之摺積法
重疊-相加算法,則是將子結果 y_L 的後段,
與下一段子結果 y_L 的前段,
像拼圖一樣根據對應的編號進行相加,
重建出完整母結果 y 的演算法。

具體來說,重疊-儲存算法切分問題的主要考量為輸入 x ,
在實作這個算法的時候必須先決定將輸入 x 切分的長度。

假設決定將輸入切分為長度 L 的區間,
我們可以直接寫出所需要的切分區段描述:

x_{Lk}[n] = x[n] (n \in kL \sim ((k+1)L - 1) , k \in 0 \sim (\left \lceil \frac{N}{L} \right \rceil - 1))

以及結果的描述:

y[n] = x[n]h[n] = \sum_{k=0}^{(\left \lceil \frac{N}{L} \right \rceil - 1)} x_{Lk}[n]h[n] = \sum_{k=0}^{(\left \lceil \frac{N}{L} \right \rceil - 1)} y_{Lk}[n]

這個算法特別需要注意的便是,重建出的母結果 y 牽涉到子結果對應項目的相加,要正確地對上才會是正確的結果,例如 y[L+1] 就會是兩個子結果的相加:
y[L+1] = \underbrace{ x[(L+1) - (M-1)]h[M+1] + \dots + x[(L+1) - 2]h[2] }_{y_L0[L+1]} + \underbrace{ x[(L+1) - 1]h[1] + x[(L+1) - 0]h[0] }_{y_L1[L+1]}

而重疊-相加算法便是得名於其在輸出的子結果上會有項目編號重疊的部分,須將其進行相加才能得到母結果。

這個算法的好處在於切分的方式相當直覺,而且在同樣的分段數量下,每個子問題需要處理的訊號長度較短,子問題計算量較小。

但缺點為可能需要一些額外的加法,而且就算採取直接計算卷積的方式來計算子問題,也仍需要補零在訊號後方。

複雜度與比較
雖然區塊卷積具有兩種差異不小的算法,但總體來說,對於每個固定架構的子問題,所需要的計算量都是常數,
而子問題的數量則正比於輸入訊號長度 N ,所以計算複雜度對兩個算法來說都是\Omicron(N)\!,
相較於使用快速傅里叶变换技巧進行整個問題的計算複雜度\Omicron(N \log N)\!要來的低。

所以一般來說在 N 很大時,顯然的使用區塊卷積會較處理整個問題來的有效率。

不過,有時候輸入訊號 x 並不具有足夠的長度可以保證區塊卷積較有效率,
又或者在一些極端情況下,直接按定義進行計算卷積可能會更有效率,
所以接著便是要稍微談論一下不同主流計算卷積方法的計算量,
以及如何根據輸入訊號長度 N 與固定訊號長度 M 間的關係選擇算法。

直接計算卷積
如果是根據卷積的定義,直接以乘法及加法計算卷積的話,
我們可以知道,每一個 x[n] 都需要輪流和每一個 h[n] 進行相乘,
所以總共需要的乘法數量為

N \times M

又因為訊號可能是複數,所以上述的乘法為複數乘法,將其統一以時數乘法進行實作後可以得到需要的實數乘法數為

3N \times M

可以發現,其實直接計算卷積對於輸入長度 N 來說,其複雜度也是\Omicron(N)\!,
但其常數為 3M ,當 M 增大時,運算量會增加很快,運算效率會不及區塊卷積,
但相對的當 M 落在很小的值時,直接計算卷積可能會是較高效率的計算方法。

使用快速傅立葉變換技巧計算卷積
除了直接計算卷積之外,另外一種主流的卷積計算方法,便是透過卷积定理,
將卷積的計算化為兩次傅立葉變換相乘後再進行反傅立葉變換的問題,
實作如下:

y[n] = x[n]*h[n] = IFFT_P(FFT_P(x[n]) \times FFT_P(h[n]))

其中FFT_P為P點快速傅立葉變換,
IFFT_P為逆變換(形式與快速傅立葉變換相同),P根據圓周卷積定理必須大於(M+N-1),
而 x[n] 及 h[n] 需透過補零來補齊長度。

如果採用這個方法,並且假設 h[n] 為固定不變,可以將 FFT_P(h[n]) 計算一次後儲存備用,
那可以得到所需乘法數為:

2 MUL_P + 3P

其中MUL_P為P點快速傅立葉變換所需的乘法數,
估計為 \frac{3}{2} P \log_2 P,
乘兩倍是因為需計算 FFT_P(x[n]) 與最後的逆變換,
3P則是因為FFT_P(x[n]) \times FFT_P(h[n]) 含有P個複數的乘法,統一轉換為實數乘法則是3P。

將 P=(M+N-1) \approx (M+N)與MUL_P的估計值帶入上面的所需乘法數進行估計,可以得到所需乘法數為:

3(M+N)(\log_2 (M+N) + 1) \approx 3(M+N)\log_2 (M+N)

可以發現,使用快速傅立葉變換技巧計算卷積對於輸入長度 N 來說,其複雜度是\Omicron(N \log N)\!,
但當 M 大約與 N 落在接近的數量級時,運算量會較其他方法低上許多。

使用區塊卷積搭配快速傅立葉算法計算卷積
最後,另外一種主流的卷積計算法便是先利用區塊卷積將問題拆分,再將拆分後的問題透過快速傅立葉算法來處理,
而之所以通常不是搭配直接計算的方式,
是因為我們通常會假設拆分過後的問題中,要進行卷積的兩訊號長度會差不多長,
而根據前面的討論,我們已經知道使用快速傅立葉轉換來進行計算會有效率的多。

而使用區塊卷積搭配快速傅立葉算法計算卷積的複雜度,
和前面兩種方法的計算複雜度只與M , N有關不同,
分段所使用的段落長度L也會很大程度的影響這個方法的計算量。

以重疊-相加算法作為代表,
假設採用的分段長度為L,
可以得到分段的段數S為:

S = \left \lceil \frac{N}{L} \right \rceil \approx \frac{N}{L}

而每個子問題的實際計算量根據上面可得:

2 MUL_P + 3 P

乘上需要解的子問題數可得實際總計算量:

2S \times MUL_P + 3S \times P

同樣根據上面可再進一步估計得:

\frac{N}{L} 3(M+L-1)(\log_2 (M+L-1) + 1) \approx \frac{N}{L} 3(M+L)\log_2 (M+L)

可以發現這個計算法的複雜度為\Omicron(N)\!,
並且常數約為\frac{3(M+L)}{L} \log_2 (M+L) ,
在M的數值略大時會較採取直接計算會較有效率,
但在M與N的數量級接近時,
因為M的大小被迫推升,
導致比起直接對整段訊號進行快速傅立葉算法的表現還要差。

所以根據上述可以估計這個方法可以展現優勢的區段,
大概是界在前面兩種方法的M , N比例之間。

而在實際應用時,因為快速傅立葉轉換的計算量在局部會有不小的波動,並不如理論值估計那般平滑,所以為了得到較好的計算效率,
除了須將上面的計算量視作L的函數,透過數值繪圖的方法或微分,求得在理論估計上最適合的L值之外,
還需在求得的最適合子問題大小附近進行找尋,
挑選數個可能的合適子問題傅立葉轉換大小,
再將反求得的分段長度帶回計算運算量,
以求得真正的最佳分段長度L。

例如:在

N = 2000,M=17

的情況下,求得的理論最適合分段長度為

L_0 = 85

理論最適合子問題傅立葉轉換大小為

P_0 = L_0 + M -1 = 101點

但在72點傅立葉轉換的計算量有一個明顯的低谷,將由此逆推得到的可能分段長度

L = 72 - M + 1 = 56

與其他可能值都帶回去計算實際計算量後,可以確認56才是最佳的問題分段長度。

最後,同樣是基於快速傅立葉轉換的計算量在局部的波動,在極為特殊的情況下,如果M \leq 4,
我們便有可能使子問題僅需使用4點傅立葉轉換,而在這樣的情況下雖會使用大量的加法,
但可不必在傅立葉轉換上使用任何乘法,可能會較直接計算卷積有效率。

參考資料
Jian-Jiun Ding, [http://djj.ee.ntu.edu.tw/ADSP.htm advanced digital signal processing lecture note] , the Department of Electrical Engineering, National Taiwan University (NTU), Taipei, Taiwan, 2022.

外部連結
重疊-儲存之摺積法

重疊-相加之摺積法

评论 (0)

  • 还没有评论,来抢沙发吧。