捲積

在泛函分析中，**捲積**（又稱**疊積**（convolution）、**褶積**或**旋積**），是透過兩個函數 _f_ 和 _g_ 生成第三個函數的一種數學算子，表徵函數 _f_ 與經過翻轉和平移的 _g_ 的乘積函數所圍成的曲邊梯形的面積。如果將參加摺積的一個函數看作區間的指示函數，摺積還可以被看作是「滑動平均」的推廣。

## 簡單介紹

摺積是數學分析中一種重要的運算。設： $f(x)$ 、 $g(x)$ 是 $\mathbb {R} $ 上的兩個可積函數，作積分：

　
: $\int _{-\infty }^{\infty }f(\tau )g(x-\tau )\,\mathrm {d} \tau $

可以證明，關於幾乎所有的 $x\in (-\infty ,\infty )$ ，上述積分是存在的。這樣，隨著 $x$ 的不同取值，這個積分就定義了一個新函數 $h(x)$ ，稱為函數 $f$ 與 $g$ 的摺積，記為 $h(x)=(f*g)(x)$ 。我們可以輕易驗證： $(f*g)(x)=(g*f)(x)$ ，並且 $(f*g)(x)$ 仍為可積函數。這就是說，把摺積代替乘法， $L^{1}(R^{1})$ 空間是一個代數，甚至是巴拿赫代數。雖然這裡為了方便我們假設 $\textstyle f,g\in L^{1}(\mathbb {R} )$ ，不過捲積只是運算符號，理論上並不需要對函數 $f,g$ 有特別的限制，雖然常常要求 $f,g$ 至少是可測函數（measurable function）（如果不是可測函數的話，積分可能根本沒有意義），至於生成的摺積函數性質會在運算之後討論。

摺積與傅立葉轉換有著密切的關係。例如兩函數的傅立葉轉換的乘積等於它們摺積後的傅立葉轉換，利用此一性質，能簡化傅立葉分析中的許多問題。

由摺積得到的函數 $f*g$ 一般要比 $f$ 和 $g$ 都光滑。特別當 $g$ 為具有緊支集的光滑函數， $f$ 為局部可積時，它們的摺積 $f*g$ 也是光滑函數。利用這一性質，對於任意的可積函數 $f$ ，都可以簡單地構造出一列逼近於 $f$ 的光滑函數列 $f_{s}$ ，這種方法稱為函數的光滑化或正則化。

摺積的概念還可以推廣到數列、測度以及廣義函數上去。

## 定義

函數 $f,g$ 是定義在 $\mathbb {R} ^{n}$ 上的可測函數（measurable function）， $f$ 與 $g$ 的摺積記作 $f*g$ ，它是其中一個函數反轉（reverse），並平移後，與另一個函數的乘積的積分，是一個對平移量的函數，也就是：

　
: $(f*g)(t)\ \ \,{\stackrel {\mathrm {def} }{=}}\ \int _{\mathbb {R} ^{n}}f(\tau )g(t-\tau )\,d\tau $ 。

如果函數不是定義在 $\mathbb {R} ^{n}$ 上，可以把函數定義域以外的值都規定成零，這樣就變成一個定義在 $\mathbb {R} ^{n}$ 上的函數。

## 離散摺積

對於定義在整數 $\mathbb {Z} $ 上的函數 $f,g$ ，摺積定義為

　
: $(f*g)[n]\ \ {\stackrel {\mathrm {def} }{=}}\ \sum _{m=-\infty }^{\infty }{f[m]g[n-m]}$

　
: : : $=\sum _{m=-\infty }^{\infty }f[n-m]\,g[m].$

這裡一樣把函數定義域以外的值當成零，所以可以擴展函數到所有整數上（如果本來不是的話）。

當 $g[n]$ 的支撐集（support）為有限長度 $M$ ，上式會變成有限和：

　
: $(f*g)[n]=\sum _{m=-M}^{M}f[n-m]g[m].$

### 計算離散摺積的方法

計算摺積 $f[n]*g[n]$ 有三種主要的方法，分別為

1. 直接計算（Direct Method）
2. 快速傅立葉轉換（FFT）
3. 分段摺積（sectioned convolution）

方法 1 是直接利用定義來計算摺積，而方法 2 和 3 都是用到了 FFT 來快速計算摺積。也有不需要用到 FFT 的作法，如使用數論轉換。

#### 方法 1：直接計算

* 作法：利用摺積的定義

　
: : $y[n]=f[n]*g[n]=\sum _{m=0}^{M-1}f[n-m]g[m]$

* 若 $f[n]$ 和 $g[n]$ 皆為實數訊號，則需要 $MN$ 個乘法。
* 若 $f[n]$ 和 $g[n]$ 皆為更一般性的複數訊號，不使用複數乘法的快速演算法，會需要 $4MN$ 個乘法；但若使用複數乘法的快速演算法，則可簡化至 $3MN$ 個乘法。

　
: 因此，使用定義直接計算摺積的複雜度為 $O(MN)$ 。

#### 方法 2：快速傅立葉轉換（FFT）

* 概念：由於兩個離散訊號在時域（time domain）做摺積相當於這兩個訊號的離散傅立葉轉換在頻域（frequency domain）做相乘：

　
: : $y[n]=f[n]*g[n]\leftrightarrow Y[f]=F[f]G[f]$

　
: ，可以看出在頻域的計算較簡單。

* 作法：因此這個方法即是先將訊號從時域轉成頻域：

　
: : $F[f]=DFT_{P}(f[n]),G[f]=DFT_{P}(g[n])$

　
: ，於是

　
: : $Y[f]=DFT_{P}(f[n])DFT_{P}(g[n])$

　
: ，最後再將頻域訊號轉回時域，就完成了摺積的計算：

　
: : $y[n]=IDFT_{P}{DFT_{P}(f[n])DFT_{P}(g[n])}$

　
: 總共做了 2 次 DFT 和 1 次 IDFT。

* 特別注意 DFT 和 IDFT 的點數 $P$ 要滿足 $P\geq M+N-1$ 。
* 由於 DFT 有快速演算法 FFT，所以運算量為 $O(P\log _{2}P)$
* 假設 $P$ 點 DFT 的乘法量為 $a$ ， $f[n]$ 和 $g[n]$ 為一般性的複數訊號，並使用複數乘法的快速演算法，則共需要 $3a+3P$ 個乘法。

#### 方法 3：分段摺積（sectioned convolution）

* 概念：將 $f[n]$ 切成好幾段，每一段分別和 $g[n]$ 做摺積後，再將結果相加。
* 作法：先將 $f[n]$ 切成每段長度為 $L$ 的區段（ $L>M$ ），假設共切成 S 段：

　
: : $f[n](n=0,1,...,N-1)\to f_{1}[n],f_{2}[n],f_{3}[n],...,f_{S}[n](S=\left\lceil {\frac {N}{L}}\right\rceil )$

　
: : Section 1: $f_{1}[n]=f[n],n=0,1,...,L-1$

　
: : Section 2: $f_{2}[n]=f[n+L],n=0,1,...,L-1$

　
: : : $\vdots $

　
: : Section r: $f_{r}[n]=f[n+(r-1)L],n=0,1,...,L-1$

　
: : : $\vdots $

　
: : Section S: $f_{S}[n]=f[n+(S-1)L],n=0,1,...,L-1$

　
: ， $f[n]$ 為各個 section 的和

　
: : $f[n]=\sum _{r=1}^{S}f_{r}[n+(r-1)L]$ 。

　
: 因此，

　
: : $y[n]=f[n]*g[n]=\sum _{r=1}^{S}\sum _{m=0}^{M-1}f_{r}[n+(r-1)L-m]g[m]$ ，

　
: 每一小段作摺積則是採用方法 2，先將時域訊號轉到頻域相乘，再轉回時域：

　
: : $y[n]=IDFT(\sum _{r=1}^{S}\sum _{m=0}^{M-1}DFT_{P}(f_{r}[n+(r-1)L-m])DFT_{P}(g[m])),P\geq M+L-1$ 。

* 總共只需要做 $P$ 點 FFT $2S+1$ 次，因為 $g[n]$ 只需要做一次 FFT。
* 假設 $P$ 點 DFT 的乘法量為 $a$ ， $f[n]$ 和 $g[n]$ 為一般性的複數訊號，並使用複數乘法的快速演算法，則共需要 $(2S+1)a+3SP$ 個乘法。
* 運算量： ${\frac {N}{L}}3(L+M-1)[\log _{2}(L+M-1)+1]$
* 運算複雜度： $O(N)$ ，和 $N$ 呈線性，較方法 2 小。
* 分為 Overlap-Add 和 Overlap-Save 兩種方法。

**分段摺積: Overlap-Add**

欲做 $x[n]*h[n]$ 的分段摺積分， $x[n]$ 長度為 $N$ ， $h[n]$ 長度為 $M$ ,

Step 1: 將 $x[n]$ 每 $L$ 分成一段

Step 2: 再每段 $L$ 點後面添加 $M-1$ 個零，變成長度 $L+M-1$

Step 3: 把 $h[n]$ 添加 $L-1$ 個零，變成長度 $L+M-1$ 的 $h'[n]$

Step 4: 把每個 $x[n]$ 的小段和 $h'[n]$ 做快速摺積，也就是 $IDFT_{L+M-1}\{{DFT_{L+M-1}(x[n])DFT_{L+M-1}(h'[n])}\}$ ，每小段會得到長度 $L+M-1$ 的時域訊號

Step 5: 放置第 $i$ 個小段的起點在位置 $L\times i$ 上， $i=0,1,...,\lceil {\frac {N}{L}}\rceil -1$

Step 6: 會發現在每一段的後面 $M-1$ 點有重疊，將所有點都相加起來，顧名思義 Overlap-Add，最後得到結果

舉例來說:

$$
x[n]=[1,2,3,4,5,-1,-2,-3,-4,-5,1,2,3,4,5]$ , 長度 $N=15
$$

$$
h[n]=[1,2,3]$ , 長度 $M=3
$$

令 $L=5$ 

* 

令 $L=5$ 切成三段，分別為 $x_{0}[n],x_{1}[n],x_{2}[n]$ , 每段填 $M-1$ 個零，並將 $h[n]$ 填零至長度 $L+M-1$

* 

將每一段做 $IDFT_{L+M-1}\{{DFT_{L+M-1}(x[n])DFT_{L+M-1}(h'[n])}\}$

* 

若將每小段擺在一起，可以注意到第一段的範圍是 $0\thicksim 6$ ，第二段的範圍是 $5\thicksim 11$ ，第三段的範圍是 $10\thicksim 16$ ，三段的範圍是有重疊的

* 

最後將三小段加在一起，並將結果和未分段的摺積做比較，上圖是分段的結果，下圖是沒有分段並利用快速摺積所算出的結果，驗證兩者運算結果相同。

* 

**分段摺積: Overlap-Save**

欲做 $x[n]*h[n]$ 的分段摺積分， $x[n]$ 長度為 $N$ ， $h[n]$ 長度為 $M$ ,

Step 1: 將 $x[n]$ 前面填 $M-1$ 個零

Step 2: 第一段 $i=0$ , 從新的 $x[n]$ 中 $L\times i-(M-1)\times i$ 取到 $L\times (i+1)-(M-1)\times i-1$ 總共 $L$ 點當做一段，因此每小段會重複取到前一小段的 $M-1$ 點，取到新的一段全為零為止

Step 3: 把 $h[n]$ 添加 $L-M$ 個零，變成長度 $L$ 的 $h'[n]$

Step 4: 把每個 $x[n]$ 的小段和 $h'[n]$ 做快速摺積，也就是 $IDFT_{L}\{{DFT_{L}(x[n])DFT_{L}(h'[n])}\}$ ，每小段會得到長度 $L$ 的時域訊號

Step 5: 對於每個 $i$ 小段，只會保留末端的 $L-(M-1)$ 點，因此得名 Overlap-Save

Step 6: 將所有保留的點合再一起，得到最後結果

舉例來說:

$$
x[n]=[1,2,3,4,5,6,7,8,9,10,11,12,13,14,15]$ , 長度 $N=15
$$

$$
h[n]=[1,2,3]$ , 長度 $M=3
$$

令 $L=7$ 

* 

將 $x[n]$ 前面填 $M-1$ 個零以後，按照 Step 2 的方式分段，可以看到每一段都重複上一段的 $M-1$ 點

* 

再將每一段做 $IDFT_{L}\{{DFT_{L}(x[n])DFT_{L}(h'[n])}\}$ 以後可以得到

* 

保留每一段末端的 $L-(M-1)$ 點，擺在一起以後，可以注意到第一段的範圍是 $0\thicksim 4$ ，第二段的範圍是 $5\thicksim 9$ ，第三段的範圍是 $10\thicksim 14$ ，第四段的範圍是 $15\thicksim 16$ ，四段的範圍是沒有重疊的

* 

將結果和未分段的摺積做比較，下圖是分段的結果，上圖是沒有分段並利用快速摺積所算出的結果，驗證兩者運算結果相同。

* 

至於為什麼要把前面 $M-1$ 丟掉？

以下以一例子來闡述:

 $x[n]=[1,2,3,4,5,6,7,8,9,10]$ , 長度 $L=10$ ，

 $h[n]=[1,2,3,4,5]$ , 長度 $M=5$ ,

第一條藍線代表 $y$ 軸，而兩條藍線之間代表長度 $L$ ，是在做快速摺積時的週期

* 

當在做快速摺積時 $IDFT_{L}\{{DFT_{L}(x[n])DFT_{L}(h'[n])}\}$ ，是把訊號視為週期 $L$ ，在時域上為循環摺積分，

而在一開始前 $M-1$ 點所得到的值，是 $h[0],h[6],h[7],h[8],h[9]$ 和 $x[0],x[6],x[7],x[8],x[9]$ 內積的值，

然而 $h[6],h[7],h[8],h[9]$ 這 $M-1$ 個值應該要為零，以往在做快速摺積時長度為 $L+M-1$ 時不會遇到這些問題，

而今天因為在做快速摺積時長度為 $L$ 才會把這 $M-1$ 點算進來，因此我們要丟棄這 $M-1$ 點內積的結果

* 

為了要丟棄這 $M-1$ 點內積的結果，位移 $h[-n]$ $M-1$ 點，並把位移以後內積合的值才算有效。

* 

#### 應用時機

以上三種方法皆可用來計算摺積，其差別在於所需母體乘法量不同。基於運算量以及效率的考量，在計算摺積時，通常會選擇所需母體乘法量較少的方法。

以下根據 $f[n]$ 和 $g[n]$ 的長度（ $N,M$ ）分成 5 類，並列出適合使用的方法：

 1. $M$ 為一非常小的整數 - 直接計算
2. $M\ll N$ - 分段摺積
3. $M\approx N$ - 快速傅立葉轉換
4. $M\gg N$ - 分段摺積
5. $N$ 為一非常小的整數 - 直接計算

基本上，以上只是粗略的分類。在實際應用時，最好還是算出三種方法所需的總乘法量，再選擇其中最有效率的方法來計算摺積。

##### 例子

Q1：當 $N=2000,M=17$ ，適合用哪種方法計算摺積？

Ans：

　
: 方法 1：所需乘法量為 $3MN=102000$

　
: 方法 2： $P\geq M+N-1=2016$ ，而 2016 點的 DFT 最少乘法數 $a=12728$ ，所以總乘法量為 $3(a+P)=44232$

　
: 方法 3：

　
: 若切成 8 塊（ $S=8$ ），則 $L=250,P\geq M+L-1=266$ 。選 $P=288$ ，則總乘法量為 $(2S+1)a+3SP=26632$ ，比方法 1 和 2 少了很多。

　
: 但是若要找到最少的乘法量，必須依照以下步驟

　
: (1) 先找出 $L$ : 解 $L$  : ${\frac {\partial {{\frac {N}{L}}3(L+M-1)[\log _{2}(L+M-1)+1]}}{\partial L}}=0$

　
: (2) 由 $P\geq L+M-1$ 算出點數在 $P$ 附近的 DFT 所需最少的乘法量，選擇 DFT 的點數

　
: (3) 最後由 $L=P+1-M$ 算出 $L_{opt}$

　
: 因此，

　
: (1) 由運算量對 $L$ 的偏微分為 0 而求出 $L=85$

　
: (2) $P\geq L+M-1=101$ ，所以選擇 101 點 DFT 附近點數乘法量最少的點數 $P=96$ 或 $P=120$ 。

　
: (3-1) 當 $P=96\to a=280,L=P+1-M=80\to S=25$ ，總乘法量為 $(2S+1)a+3SP=21480$ 。

　
: (3-2) 當 $P=120\to a=380,L=P+1-M=104\to S=20$ ，總乘法量為 $(2S+1)a+3SP=22780$ 。

　
: 由此可知，切成 20 塊會有較好的效率，而所需總乘法量為 21480。

* 因此，當 $N=2000,M=17$ ，所需總乘法量：分段摺積 < 快速傅立葉轉換 < 直接計算。故，此時選擇使用分段摺積來計算摺積最適合。

Q2：當 $N=1024,M=3$ ，適合用哪種方法計算摺積？

Ans：

　
: 方法 1：所需乘法量為 $3MN=9216$

　
: 方法 2： $P\geq M+N-1=1026$ ，選擇 1026 點 DFT 附近點數乘法量最少的點數， $\to P=1152,a=7088$ 。

　
: : 因此，所需乘法量為 $3(a+P)=24342$

　
: 方法 3：

　
: : (1) 由運算量對 $L$ 的偏微分為 0 而求出 $L=5$

　
: (2) $P\geq L+M-1=7$ ，所以選擇 7 點 DFT 附近點數乘法量最少的點數 $P=8$ 或 $P=6$ 或 $P=4$ 。

　
: (3-1) 當 $P=8\to a=4,L=P+1-M=6\to S=171$ ，總乘法量為 $(2S+1)a+3SP=5476$ 。

　
: (3-2) 當 $P=6\to a=4,L=P+1-M=4\to S=256$ ，總乘法量為 $(2S+1)a+3SP=6660$ 。

　
: (3-3) 當 $P=4\to a=0,L=P+1-M=2\to S=512$ ，總乘法量為 $(2S+1)a+3SP=6144$ 。

　
: 由此可知，切成 171 塊會有較好的效率，而所需總乘法量為 5476。

* 因此，當 $N=1024,M=3$ ，所需總乘法量：分段摺積 < 直接計算 < 快速傅立葉轉換。故，此時選擇使用分段摺積來計算摺積最適合。
* 雖然當 $M$ 是個很小的正整數時，大致上適合使用直接計算。但實際上還是將 3 個方法所需的乘法量都算出來，才能知道用哪種方法可以達到最高的效率。

Q3：當 $N=1024,M=600$ ，適合用哪種方法計算摺積？

Ans：

　
: 方法 1：所需乘法量為 $3MN=1843200$

　
: 方法 2： $P\geq M+N-1=1623$ ，選擇 1026 點 DFT 附近點數乘法量最少的點數， $\to P=2016,a=12728$ 。

　
: : 因此，所需乘法量為 $3(a+P)=44232$

　
: 方法 3：

　
: : (1) 由運算量對 $L$ 的偏微分為 0 而求出 $L=1024$

　
: (2) $P\geq L+M-1=1623$ ，所以選擇 1623 點 DFT 附近點數乘法量最少的點數 $P=2016$ 。

　
: (3) 當 $P=2016\to a=12728,L=P+1-M=1417\to S=1$ ，總乘法量為 $(2S+1)a+3SP=44232$ 。

　
: 由此可知，此時切成一段，就跟方法 2 一樣，所需總乘法量為 44232。

* 因此，當 $N=1024,M=600$ ，所需總乘法量：快速傅立葉轉換 = 分段摺積 < 直接計算。故，此時選擇使用分段摺積來計算摺積最適合。

## 多元函數摺積

按照翻轉、平移、積分的定義，還可以類似的定義多元函數上的積分：

　
: $(f*g)(t_{1},t_{2},\cdots ,t_{n})=\int \int \cdots \int f(\tau _{1},\tau _{2},\cdots ,\tau _{n})g(t_{1}-\tau _{1},t_{2}-\tau _{2},\cdots ,t_{n}-\tau _{n},)\,d\tau _{1}d\tau _{2}\cdots d\tau _{n}$

## 性質

各種摺積算子都滿足下列性質：

**交換律**
: $f*g=g*f\,$

**結合律**
: $f*(g*h)=(f*g)*h\,$

**分配律**
: $f*(g+h)=(f*g)+(f*h)\,$

**數乘結合律**
: $a(f*g)=(af)*g=f*(ag)\,$

其中 $a$ 為任意實數（或複數）。

**微分定理**
: ${\mathcal {D}}(f*g)={\mathcal {D}}f*g=f*{\mathcal {D}}g\,$

其中 D_f_ 表示 _f_ 的微分，如果在離散域中則是指差分算子，包括前向差分與後向差分兩種：

* 前向差分： ${\mathcal {D}}^{+}f(n)=f(n+1)-f(n)$
* 後向差分： ${\mathcal {D}}^{-}f(n)=f(n)-f(n-1)$

## 摺積定理

**摺積定理**指出，函數摺積的傅立葉轉換是函數傅立葉轉換的乘積。即，一個域中的摺積相當於另一個域中的乘積，例如時域中的摺積就對應於頻域中的乘積。

　
: ${\mathcal {F}}(f*g)={\mathcal {F}}(f)\cdot {\mathcal {F}}(g)$

其中 ${\mathcal {F}}(f)$ 表示 _f_ 的傅立葉轉換。

這一定理對拉普拉斯轉換、雙邊拉普拉斯轉換、Z 轉換、Mellin 轉換和 Hartley 轉換（參見 Mellin inversion theorem）等各種傅立葉轉換的變體同樣成立。在調和分析中還可以推廣到在局部緊緻的阿貝爾群上定義的傅立葉轉換。

利用摺積定理可以簡化摺積的運算量。對於長度為 $n$ 的序列，按照摺積的定義進行計算，需要做 $2n-1$ 組對位乘法，其計算複雜度為 ${\mathcal {O}}(n^{2})$ ；而利用傅立葉轉換將序列轉換到頻域上後，只需要一組對位乘法，利用傅立葉轉換的快速算法之後，總的計算複雜度為 ${\mathcal {O}}(n\log n)$ 。這一結果可以在快速乘法計算中得到應用。

## 在群上的摺積

若 _G_ 是有某 _m_ 測度的群（例如郝斯多夫空間上哈爾測度下局部緊緻的拓撲群），對於 _G_ 上 _m_- 勒貝格可積的實數或複數函數 _f_ 和 _g_，可定義它們的摺積：

　
: $(f*g)(x)=\int _{G}f(y)g(xy^{-1})\,dm(y)\,$

對於這些群上定義的摺積同樣可以給出諸如摺積定理等性質，但是這需要對這些群的表示理論以及調和分析的彼得 - 外爾定理。

## 應用

摺積在科學、工程和數學上都有很多應用：

* 代數中，整數乘法和多項式乘法都是摺積。
* 圖像處理中，用作圖像模糊、銳化、邊際檢測。
* 統計學中，加權的滑動平均是一種摺積。
* 機率論中，兩個統計獨立變數 X 與 Y 的和的機率密度函數是 X 與 Y 的機率密度函數的摺積。
* 聲學中，回聲可以用源聲與一個反映各種反射效應的函數的摺積表示。
* 電子工程與訊號處理中，任一個線性系統的輸出都可以通過將輸入訊號與系統函數（系統的脈衝響應）做摺積獲得。
* 物理學中，任何一個線性系統（符合疊加原理）都存在摺積。

## 參見

* 反摺積
* 自我相關函數
* 傅立葉轉換

## 外部連結

* Convolution. PlanetMath.
* Visual convolution Java Applet （頁面存檔備份，存於網際網路檔案館）
* Lecture notes, Jian-Jiun Ding (2013), Advanced Digital Signal Processing（頁面存檔備份，存於網際網路檔案館）