图像分割--边缘检测算子
一阶微分算子
边缘可以定位为一阶导数较大的像素处
Robert 算子
第一个边缘检测算子,1963年 Lawrence Roberts 提出
算子模板及差分公式
(a)
G
x
=
△
x
f
(
x
,
y
)
=
f
(
x
+
1
,
y
+
1
)
−
f
(
x
,
y
)
G_x=\bigtriangleup_xf(x,y)=f(x+1,y+1)-f(x,y)
Gx=△xf(x,y)=f(x+1,y+1)−f(x,y)
(b)
G
y
=
△
y
f
(
x
,
y
)
=
f
(
x
+
1
,
y
)
−
f
(
x
,
y
+
1
)
G_y=\bigtriangleup_yf(x,y)=f(x+1,y)-f(x,y+1)
Gy=△yf(x,y)=f(x+1,y)−f(x,y+1)
Sobel算子
1968年提出,1973年出现在专著《Pattern Classification and Scene Analysis》中
算子模板及差分公式
(a)
G
x
=
△
x
f
(
x
,
y
)
=
f
(
x
−
1
,
y
−
1
)
+
2
f
(
x
−
1
,
y
)
+
f
(
x
−
1
,
y
−
1
)
−
[
f
(
x
+
1
,
y
+
1
)
+
2
f
(
x
+
1
,
y
)
+
f
(
x
+
1
,
y
+
1
)
]
G_x=\bigtriangleup_xf(x,y)=f(x-1,y-1)+2f(x-1,y)+f(x-1,y-1)-[f(x+1,y+1)+2f(x+1,y)+f(x+1,y+1)]
Gx=△xf(x,y)=f(x−1,y−1)+2f(x−1,y)+f(x−1,y−1)−[f(x+1,y+1)+2f(x+1,y)+f(x+1,y+1)]
(b) G y = △ y f ( x , y ) = f ( x + 1 , y + 1 ) + 2 f ( x , y + 1 ) + f ( x − 1 , y + 1 ) − [ f ( x + 1 , y − 1 ) + 2 f ( x , y − 1 ) + f ( x − 1 , y − 1 ) ] G_y=\bigtriangleup_yf(x,y)=f(x+1,y+1)+2f(x,y+1)+f(x-1,y+1)-[f(x+1,y-1)+2f(x,y-1)+f(x-1,y-1)] Gy=△yf(x,y)=f(x+1,y+1)+2f(x,y+1)+f(x−1,y+1)−[f(x+1,y−1)+2f(x,y−1)+f(x−1,y−1)]
Sobel算子的变种–各向同性Sobel算子
各向同性 Sobel 算子权值更准确,按 距离越远影响越小 来说。
Prewitt 算子
1970 年,《Object Enhancement and Extraction》in Picture processing and Psychopictorics
(a)
G
x
=
△
x
f
(
x
,
y
)
=
f
(
x
−
1
,
y
−
1
)
+
f
(
x
−
1
,
y
)
+
f
(
x
−
1
,
y
−
1
)
−
[
f
(
x
+
1
,
y
+
1
)
+
f
(
x
+
1
,
y
)
+
f
(
x
+
1
,
y
+
1
)
]
G_x=\bigtriangleup_xf(x,y)=f(x-1,y-1)+f(x-1,y)+f(x-1,y-1)-[f(x+1,y+1)+f(x+1,y)+f(x+1,y+1)]
Gx=△xf(x,y)=f(x−1,y−1)+f(x−1,y)+f(x−1,y−1)−[f(x+1,y+1)+f(x+1,y)+f(x+1,y+1)]
(b) G y = △ y f ( x , y ) = f ( x + 1 , y + 1 ) + f ( x , y + 1 ) + f ( x − 1 , y + 1 ) − [ f ( x + 1 , y − 1 ) + f ( x , y − 1 ) + f ( x − 1 , y − 1 ) ] G_y=\bigtriangleup_yf(x,y)=f(x+1,y+1)+f(x,y+1)+f(x-1,y+1)-[f(x+1,y-1)+f(x,y-1)+f(x-1,y-1)] Gy=△yf(x,y)=f(x+1,y+1)+f(x,y+1)+f(x−1,y+1)−[f(x+1,y−1)+f(x,y−1)+f(x−1,y−1)]
算子使用方式及实现
算子模板在图像上滑动,覆盖图像做卷积后得到检测幅值
幅值:
G
=
G
x
2
+
G
y
2
G = \sqrt{G_x^2+G_y^2}
G=Gx2+Gy2 或
G
=
m
a
x
(
∣
G
x
∣
,
∣
G
y
∣
)
G=max(|G_x|, |G_y|)
G=max(∣Gx∣,∣Gy∣) 或
G
=
∣
G
x
∣
+
∣
G
y
∣
G=|G_x|+|G_y|
G=∣Gx∣+∣Gy∣
from PIL import Image
import numpy as np
import cv2
from matplotlib import pyplot as plt
filename = '002646.jpg'
img = Image.open(filename).convert('RGB')
img_gray = np.array(img.convert('L'))
out_img = np.zeros_like(img_gray)
def robert(img_gray, out_img):
h, w = img_gray.shape
for i in range(h-1):
for j in range(w-1):
dx = img_gray[i+1][j+1]-img_gray[i][j]
dy = img_gray[i+1][j] - img_gray[i][j+1]
out_img[i][j] = np.sqrt(dx**2 + dy**2)
return out_img
def sobel(img_gray, out_img):
h, w = img_gray.shape
for i in range(h-1):
for j in range(w-1):
dx = img_gray[i-1][j+1] + 2*img_gray[i][j+1] + img_gray[i+1][j+1]-img_gray[i-1][j-1] - 2*img_gray[i][j-1] - img_gray[i+1][j-1]
dy = img_gray[i-1][j-1] + 2*img_gray[i-1][j] + img_gray[i-1][j+1] - img_gray[i+1][j-1] -2*img_gray[i+1][j] - img_gray[i+1][j+1]
out_img[i][j] = np.sqrt(dx**2 + dy**2)
return out_img
def perwitt(img_gray, out_img):
h, w = img_gray.shape
for i in range(h-1):
for j in range(w-1):
dx = img_gray[i-1][j-1]+img_gray[i][j-1]+img_gray[i+1][j-1]-img_gray[i-1][j+1]-img_gray[i][j+1]-img_gray[i+1][j+1]
dy = img_gray[i-1][j-1]+img_gray[i-1][j]+img_gray[i-1][j+1] - img_gray[i+1][j-1]-img_gray[i+1][j]-img_gray[i+1][j+1]
out_img[i][j] = np.sqrt(dx**2 + dy**2)
return out_img
out_img3 = perwitt(img_gray, out_img)
out_img = robert(img_gray, out_img)
x = cv2.Sobel(img_gray, cv2.CV_16S, 1, 0)
y = cv2.Sobel(img_gray, cv2.CV_16S, 0, 1)
absX = cv2.convertScaleAbs(x)
absY = cv2.convertScaleAbs(y)
Sobel = cv2.addWeighted(absX, 0.5, absY, 0.5, 0) ## 将两张大小相同的图片进行融合的
out_img2 = sobel(img_gray, out_img)
plt.subplot(151)
plt.imshow(img)
plt.subplot(152)
plt.imshow(out_img, cmap='gray')
plt.subplot(153)
plt.imshow(out_img2, cmap='gray')
plt.subplot(154)
plt.imshow(Sobel, cmap='gray')
plt.subplot(155)
plt.imshow(out_img3, cmap='gray')
plt.show()
## 卷积操作实现Sobel算子
sobelx = np.array([[1, 0, -1], [1, 0, -1], [1, 0, -1]])
sobely = np.array([[-1, -1, -1], [0, 0, 0], [1, 1, 1]])
x = cv2.filter2D(img_gray, cv2.CV_16S, sobelx)
y = cv2.filter2D(img_gray, cv2.CV_16S, sobely)
absX = cv2.convertScaleAbs(x)
absY = cv2.convertScaleAbs(y)
conv_sobel = cv2.addWeighted(absX, 0.5, absY, 0.5, 0)
plt.imshow(conv_sobel, cmap='gray')
plt.show()
效果:从左到右依次为原图、Robert算子、Sobel算子(差分)、Sobel算子(CV2) 、以及Perwitt算子
以及卷积操作实现的Sobel算子效果
其他的算子差分计算是一样的,这里不再附上代码。
二阶微分算子
二阶导为0的地方就是边缘
Laplace算子
是一个各项同性算子,具有旋转不变性,只关心边缘的位置而不考虑其周围的像素灰度差时比较合适,但对噪声非常敏感,通常的分割算法会将 Laplace算子与平滑算子相结合。
噪声敏感的原因
图像中噪声是高频信息,一阶导后会放大突变值,频率依旧很高,二阶导再次放大突变值,使得突变更大。对图像本身来说就是更大的噪声,故 Laplace 算子对噪声敏感。按这个推论,再加一阶导噪声更敏感。
算子模板及差分公式
上面两个模板分别为八邻域和四邻域模板
▽
2
f
(
x
,
y
)
=
[
f
(
x
+
1
,
y
)
−
f
(
x
,
y
)
]
−
[
f
(
x
,
y
)
−
f
(
x
−
1
,
y
)
]
+
[
f
(
x
,
y
+
1
)
−
f
(
x
,
y
)
]
−
[
f
(
x
,
y
)
−
f
(
x
,
y
−
1
)
]
\bigtriangledown^2f(x,y)=[f(x+1,y)-f(x,y)]-[f(x,y)-f(x-1,y)]+[f(x,y+1)-f(x,y)]-[f(x,y)-f(x,y-1)]
▽2f(x,y)=[f(x+1,y)−f(x,y)]−[f(x,y)−f(x−1,y)]+[f(x,y+1)−f(x,y)]−[f(x,y)−f(x,y−1)]
=
f
(
x
+
1
,
y
)
+
f
(
x
−
1
,
y
)
+
f
(
x
,
y
+
1
)
+
f
(
x
,
y
−
1
)
−
4
f
(
x
,
y
)
=f(x+1,y)+f(x-1,y)+f(x,y+1)+f(x,y-1)-4f(x,y)
=f(x+1,y)+f(x−1,y)+f(x,y+1)+f(x,y−1)−4f(x,y)
高斯型 Laplace 算子(Laplacian of Gaussian, LoG)
二维高斯卷积函数
G
0
(
x
,
y
)
=
1
2
π
σ
2
e
−
x
2
+
y
2
2
σ
2
G_0(x,y)=\frac{1}{2\pi \sigma ^2}e^{-\frac{x^2+y^2}{2\sigma ^2}}
G0(x,y)=2πσ21e−2σ2x2+y2
可以用这个函数,以模板中间像素为中心求出模板其他位置处的值。
具体表示
d
G
0
d
x
=
−
x
σ
2
e
−
x
2
+
y
2
2
σ
2
1
2
π
σ
2
\frac{dG_0}{dx}=-\frac{x}{\sigma ^2}e^{-\frac{x^2+y^2}{2\sigma ^2}} \frac{1}{2 \pi \sigma ^2}
dxdG0=−σ2xe−2σ2x2+y22πσ21
d
2
G
0
d
x
2
=
(
x
2
σ
4
e
−
x
2
+
y
2
2
σ
2
−
1
σ
2
e
−
x
2
+
y
2
2
σ
2
)
1
2
π
σ
2
\frac{d^2G_0}{dx^2}=(\frac{x^2}{\sigma ^4}e^{-\frac{x^2+y^2}{2\sigma ^2}}-\frac{1}{\sigma ^2}e^{-\frac{x^2+y^2}{2\sigma ^2}})\frac{1}{2 \pi \sigma ^2}
dx2d2G0=(σ4x2e−2σ2x2+y2−σ21e−2σ2x2+y2)2πσ21
△
G
0
=
(
x
2
+
y
2
σ
4
e
−
x
2
+
y
2
2
σ
2
−
1
σ
2
e
−
x
2
+
y
2
2
σ
2
)
1
2
π
σ
2
=
x
2
+
y
2
−
2
σ
2
2
π
σ
6
e
−
x
2
+
y
2
2
σ
2
\bigtriangleup G_0=(\frac{x^2+y^2}{\sigma ^4}e^{-\frac{x^2+y^2}{2\sigma ^2}}-\frac{1}{\sigma ^2}e^{-\frac{x^2+y^2}{2\sigma ^2}})\frac{1}{2 \pi \sigma ^2} \\ =\frac{x^2+y^2-2\sigma ^2}{2\pi \sigma ^6}e^{-\frac{x^2+y^2}{2\sigma ^2}}
△G0=(σ4x2+y2e−2σ2x2+y2−σ21e−2σ2x2+y2)2πσ21=2πσ6x2+y2−2σ2e−2σ2x2+y2
算子模板
高斯差分(Difference of Gaussian, DoG)
D
o
G
=
G
(
x
,
y
,
σ
1
)
−
G
(
x
,
y
,
σ
2
)
DoG=G(x,y,\sigma _1)-G(x,y, \sigma _2)
DoG=G(x,y,σ1)−G(x,y,σ2)
∂
G
∂
σ
≈
G
(
x
,
y
,
k
σ
−
G
(
x
,
y
,
σ
k
σ
−
σ
\frac{\partial G}{\partial \sigma}\approx\frac{G(x,y,k\sigma-G(x,y,\sigma}{k\sigma-\sigma}
∂σ∂G≈kσ−σG(x,y,kσ−G(x,y,σ
D
o
G
≈
(
k
−
1
)
σ
∂
G
∂
σ
DoG\approx (k-1)\sigma \frac{\partial G}{\partial \sigma}
DoG≈(k−1)σ∂σ∂G
又
∂
G
∂
σ
=
−
1
(
2
π
σ
2
)
2
4
π
σ
e
−
x
2
+
y
2
2
σ
2
+
4
π
2
π
σ
2
x
2
+
y
2
(
2
σ
2
)
2
e
−
x
2
+
y
2
2
σ
2
=
x
2
+
y
2
−
2
σ
2
2
π
σ
5
e
−
x
2
+
y
2
2
σ
2
=
σ
⋅
△
G
\frac{\partial G}{\partial \sigma}=-\frac{1}{(2\pi\sigma ^2)^2}4\pi \sigma e^{-\frac{x^2+y^2}{2\sigma ^2}}+\frac{4\pi}{2\pi \sigma^2}\frac{x^2+y^2}{(2\sigma^2)^2}e^{-\frac{x^2+y^2}{2\sigma ^2}}\\ =\frac{x^2+y^2-2\sigma ^2}{2\pi \sigma^5}e^{-\frac{x^2+y^2}{2\sigma ^2}}\\ =\sigma \cdot \bigtriangleup G
∂σ∂G=−(2πσ2)214πσe−2σ2x2+y2+2πσ24π(2σ2)2x2+y2e−2σ2x2+y2=2πσ5x2+y2−2σ2e−2σ2x2+y2=σ⋅△G
故
D
o
G
≈
(
k
−
1
)
σ
2
△
G
DoG \approx (k-1)\sigma ^2 \bigtriangleup G
DoG≈(k−1)σ2△G
即DoG与
△
G
\bigtriangleup G
△G 波形相似,但值的大小不同,故不影响边缘检测。而且从计算可以看出来,DoG计算更简单。
非微分边缘检测算子 Canny 算子
Canny 算子准确来说是一种多级边缘检测算法,具体步骤如下:
- 图像灰度化
- 高斯模糊
- 计算梯度幅值及方向
- 非最大信号压制
- 双阈值边缘连接处理
- 二值化图像输出结果
灰度化
公式
g
r
a
y
=
0.299
∗
R
+
0.587
∗
G
+
0.114
B
gray = 0.299*R+0.587*G+0.114B
gray=0.299∗R+0.587∗G+0.114B
计算梯度
可以用 Sobel 算子、Roberts 算子和 Prewitt算子等
非最大值信号抑制
用其他检测算子算出来的边缘太粗了。为了细化边界,应用非最大抑制–它只保留梯度最大位置,一般用八邻域。
具体做法
沿边界上每个像素处的梯度方向,会在它的八邻域内有两个交点。有时不在八邻域的正中点,可以线性计算交点处梯度值,即该点两侧的梯度线性组合,系数 w 是梯度方向角度的 正切值。
双阈值边缘连接处理
认为给定阈值上界和下界。像素点大于上界的认为必然是边界(强边界)。小于下界的必然不是边界。两者之间的是候选项(弱边界)。
如可以保留一半的弱边界作为边界。
滞后边界跟踪:对弱边缘的二次判断。
若弱边缘点周围八邻域存在强边缘点,则该弱边缘点变成强边缘点。
[1]Canny算子详解–代码清晰
[2] Canny边缘检测
[3]边缘算子总结
[4] Opencv-边缘检测算子及实现