索貝爾算子

**索伯算子**（Sobel operator）是圖像處理中的算子之一，有時又稱為**索伯 - 費爾德曼算子**或**索伯濾波器**，在影像處理及電腦視覺領域中常被用來做邊緣檢測。索伯算子最早是由美國計算機科學家艾爾文・索伯及蓋瑞・費爾德曼（Gary Feldman）於 1968 年在史丹佛大學的人工智慧實驗室（SAIL）所提出，因此為了表揚他們的貢獻，而用他們的名字命名。在技術上，它是一離散性差分算子，用來運算圖像亮度函數的梯度之近似值。在圖像的任何一點使用此算子，索伯算子的運算將會產生對應的梯度向量或是其範數。概念上，索伯算子就是一個小且是整數的濾波器對整張影像在水平及垂直方向上做捲積，因此它所需的運算資源相對較少，另一方面，對於影像中的頻率變化較高的地方，它所得的梯度之近似值也比較粗糙。

## 公式

該算子包含兩組 3x3 的矩陣，分別為橫向及縱向，將之與圖像作平面卷積，即按所示方程式計算：

$$
{\begin{aligned}&\left({\begin{bmatrix}a&b&c\\d&e&f\\g&h&i\end{bmatrix}}*{\begin{bmatrix}1&2&3\\4&5&6\\7&8&9\end{bmatrix}}\right)[2,2]=\\&(i\cdot 1)+(h\cdot 2)+(g\cdot 3)+\\&(f\cdot 4)+(e\cdot 5)+(d\cdot 6)+\\&(c\cdot 7)+(b\cdot 8)+(a\cdot 9)\end{aligned}}
$$

即可分別得出橫向及縱向的亮度差分近似值。如果以 $\mathbf {A} $ 代表原始圖像， $\mathbf {G_{x}} $ 及 $\mathbf {G_{y}} $ 分別代表經橫向及縱向邊緣檢測的圖像，其公式如下：

　
: $\mathbf {G_{x}} ={\begin{bmatrix}+1&0&-1\\+2&0&-2\\+1&0&-1\end{bmatrix}}*\mathbf {A} \quad {\mbox{and}}\quad \mathbf {G_{y}} ={\begin{bmatrix}+1&+2&+1\\0&0&0\\-1&-2&-1\end{bmatrix}}*\mathbf {A} $

圖像的每一個像素的橫向及縱向梯度近似值可用以下的公式結合，來計算梯度的大小。

　
: $\mathbf {G} ={\sqrt {\mathbf {G_{x}} ^{2}+\mathbf {G_{y}} ^{2}}}$

然後可用以下公式計算梯度方向。

　
: $\mathbf {\Theta } =\operatorname {arctan} \left({\mathbf {G_{y}}  \over \mathbf {G_{x}} }\right)$

以縱向邊緣為例，如果角度 $\Theta $ 等於零，即代表圖像的縱向邊緣右方較亮；若為 $\pi $ ，則左方較亮。

## 更正式來說

雖然索伯算子代表著相對沒那麼精準的影像梯度近似，但其品質也足夠在很多應用派上用場。更精確來說，對於每個點，它只用其亮度的值在周圍 3x3 的區域去近似該點的梯度，而且在近似的過程，它也只用整數去代表影像亮度的權重。

## 索伯算子在其他維度的延伸

索伯算子包含著兩個分開的運算:

* 利用三角形函數濾波器來達到平滑化 : $h(-1)=1,h(0)=2,h(1)=1$
* 一次微分的中央差值 : $h'(-1)=1,h'(0)=0,h'(1)=-1$

索伯算子在不同維度的表達形式， $x,y,z,t\in \left\{0,-1,1\right\}$ :

1D: $h_{x}'(x)=h'(x);$

2D: $h_{x}'(x,y)=h'(x)h(y)$

3D: $h_{x}'(x,y,z)=h'(x)h(y)h(z)$

4D: $h_{x}'(x,y,z,t)=h'(x)h(y)h(z)h(t)$

下面是一個三維的索伯算子在 z 軸的實例:

　
: $h_{z}'(:,:,-1)={\begin{bmatrix}+1&+2&+1\\+2&+4&+2\\+1&+2&+1\end{bmatrix}}\quad h_{z}'(:,:,0)={\begin{bmatrix}0&0&0\\0&0&0\\0&0&0\end{bmatrix}}\quad h_{z}'(:,:,1)={\begin{bmatrix}-1&-2&-1\\-2&-4&-2\\-1&-2&-1\end{bmatrix}}$

## 技術細節

索伯算子在實作上可以用很簡單的硬體和軟體：只有該點的周圍八個點需要被計算，同時在梯度向量的近似上也都只牽涉到整數的計算，而且上述的兩個離散濾波器是可分解的:

　
: ${\begin{bmatrix}1&0&-1\\2&0&-2\\1&0&-1\end{bmatrix}}={\begin{bmatrix}1\\2\\1\end{bmatrix}}{\begin{bmatrix}1&0&-1\end{bmatrix}}={\begin{bmatrix}1\\1\end{bmatrix}}*{\begin{bmatrix}1\\1\end{bmatrix}}{\begin{bmatrix}1&-1\end{bmatrix}}*{\begin{bmatrix}1&1\end{bmatrix}}$

　
: ${\begin{bmatrix}\ \ 1&\ \ 2&\ \ 1\\\ \ 0&\ \ 0&\ \ 0\\-1&-2&-1\end{bmatrix}}={\begin{bmatrix}\ \ 1\\\ \ 0\\-1\end{bmatrix}}{\begin{bmatrix}1&2&1\end{bmatrix}}={\begin{bmatrix}1\\1\end{bmatrix}}*{\begin{bmatrix}\ \ 1\\-1\end{bmatrix}}{\begin{bmatrix}1&1\end{bmatrix}}*{\begin{bmatrix}1&1\end{bmatrix}}$

因此 **G**x 和 **G**y 可推導成

　
: $\mathbf {G} _{x}={\begin{bmatrix}1\\2\\1\end{bmatrix}}*\left({\begin{bmatrix}1&0&-1\end{bmatrix}}*\mathbf {A} \right)\quad {\mbox{and}}\quad \mathbf {G} _{y}={\begin{bmatrix}\ \ 1\\\ \ 0\\-1\end{bmatrix}}*\left({\begin{bmatrix}1&2&1\end{bmatrix}}*\mathbf {A} \right)$

在特別的實作上，可分解的運算可能會有優勢因為在每個點的運算上會有少一些的算術運算。

應用捲積核 K 在像素群 P 上的虛擬碼可以表示成:

N(x,y) = Sum of { K(i,j).P(x-i,y-j)}, for i,j running from -1 to 1.

P 代表著像素矩陣，所以 N (x,y) 代表的應用捲積核 K 在像素群 P 後產生的矩陣

## 實例

經過索伯算子所產生的結果是一張二維的圖，代表的每個點的梯度。下圖中高梯度的地方會以白色的線表示。為了說明這一點，下方展示了在一張簡單的圖片上使用索伯算子的結果。

下方的灰階影像說明了梯度的方向。紅色和黃色代表梯度角度是正的，相對地，藍色和青色代表梯度角度是負的。在圓的最左邊和最右邊，垂直的梯度角度為零是因為沒有區域的變化，而在圓的最頂部和底部，水平的梯度角度為零是因為沒有區域的變化。在圓的最頂部和底部，水平梯度角度變化分別是 −π/2 和 π/2 是因為沒有區域的變化。頂部的負梯度角度代表這條邊是從亮到暗的區域界線，而底部的正梯度角度代表這條邊是從暗到亮的區域界線。其他顯示為黑色的像素是因為附近完全沒有任何變化。值得注意的是，因為角度是對應到像素值的函數，因此可能一些微小的變化也可能會導致很大的角度變化。因此，些微的雜訊可能會對角度響應變化，而這通常是不樂見的情況。由上述可知，使用角度梯度在影像處理的應用時，需要在影像去噪多加著墨以減少錯誤。

## 替代的運算子

索伯算子在轉動上沒有完美的對稱。因此沙爾想要改進這個特性。它曾提出了更大的 5x5 核心，不過後來最常被使用的是:

　
: ${\begin{bmatrix}+3&0&-3\\+10&0&-10\\+3&0&-3\end{bmatrix}}\ \ \ \ \ \ \ \ \ {\begin{bmatrix}+3&+10&+3\\0&0&0\\-3&-10&-3\end{bmatrix}}$

沙爾算子的概念源自於試圖在頻域上最小化加權的均方角度差，這個最佳化問題需要考量濾波器在數值上的一致性，因此它們是微分核。

另一個相似的改進方法是被法立德和西蒙切利所提出，它們的研究是針對更高次方的微分，相較於沙爾算子，這些濾波器沒有考量到數值的一致性。

除此之外，克龍也提出了利用數值方法的最佳化來設計微分濾波器。

2014 年，哈斯特也利用任意的三次樣條曲線去設計微分濾波器，他的研究指出可以利用長度為 7 的濾波器作雙重濾波在三次樣條曲線及三角樣條上已取得精準的一次微分和二次微分。

另外一個相似的算子是卡雅利算子，概念也是來自於索伯算子，但它是個具有完美的旋轉對稱的 3x3 捲積濾波器。

在那之後的工作，沙爾又針對了光流的預測上提出了更好的濾波器，同一年，他也提出了更好的二次微分濾波器應用在動作預測上。結果顯示使用的核心越大，就越有可能去逼近高斯濾波的微分。

## 結果比較

下圖是四張不同的梯度算子所預測出的測試圖的梯度強度

## 偽代碼

```
function sobel(A : as two dimensional image array)
	Gx=[-1 0 1; -2 0 2; -1 0 1]
	Gy=[-1 -2 -1; 0 0 0; 1 2 1]
	
	rows = size(A,1)
	columns = size(A,2)
	mag=zeros(A)

	for i=1:rows-2
		for j=1:columns-2
			S1=sum(sum(Gx.*A(i:i+2,j:j+2)))
			S2=sum(sum(Gy.*A(i:i+2,j:j+2)))

			mag(i+1,j+1)=sqrt(S1.^2+S2.^2)
		end for
	end for
	
	threshold = 70 %varies for application [0 255]
	output_image = max(mag,threshold)
	output_image(output_image==round(threshold))=0;
	return output_image
end function
```

## 參考文獻

## 參見

* 數字圖像處理
* 邊緣檢測
* 圖像梯度