Otsu’s Thresholding的工作原理
在双峰图像(像素的直方图只包含两个峰值),一个好的阈值将在这两个值的中间。Otsu算法就是从图像的像素直方图中确定了一个最佳的全局阈值
t
t
t,可以最小化weighted within-class variance
加权类内方差。
σ
w
2
(
t
)
=
q
1
(
t
)
σ
1
2
(
t
)
+
q
2
(
t
)
σ
2
2
(
t
)
\sigma_w^2(t) = q_1(t)\sigma_1^2(t)+q_2(t)\sigma_2^2(t)
σw2(t)=q1(t)σ12(t)+q2(t)σ22(t)
权重:
q
1
(
t
)
=
∑
i
=
1
t
P
(
i
)
&
q
2
(
t
)
=
∑
i
=
t
+
1
I
P
(
i
)
q_1(t) = \sum_{i=1}^{t} P(i) \quad \& \quad q_2(t) = \sum_{i=t+1}^{I} P(i)
q1(t)=i=1∑tP(i)&q2(t)=i=t+1∑IP(i)
均值:
μ
1
(
t
)
=
∑
i
=
1
t
i
P
(
i
)
q
1
(
t
)
&
μ
2
(
t
)
=
∑
i
=
t
+
1
I
i
P
(
i
)
q
2
(
t
)
\mu_1(t) = \sum_{i=1}^{t} \frac{iP(i)}{q_1(t)} \quad \& \quad \mu_2(t) = \sum_{i= t+1}^{I} \frac{iP(i)}{q_2(t)}
μ1(t)=i=1∑tq1(t)iP(i)&μ2(t)=i=t+1∑Iq2(t)iP(i)
方差:
σ
1
2
(
t
)
=
∑
i
=
1
t
[
i
−
μ
1
(
t
)
]
2
P
(
i
)
q
1
(
t
)
&
σ
2
2
(
t
)
=
∑
i
=
t
+
1
I
[
i
−
μ
2
(
t
)
]
2
P
(
i
)
q
2
(
t
)
\sigma_1^2(t) = \sum_{i=1}^{t} [i-\mu_1(t)]^2 \frac{P(i)}{q_1(t)} \quad \& \quad \sigma_2^2(t) = \sum_{i=t+1}^{I} [i-\mu_2(t)]^2 \frac{P(i)}{q_2(t)}
σ12(t)=i=1∑t[i−μ1(t)]2q1(t)P(i)&σ22(t)=i=t+1∑I[i−μ2(t)]2q2(t)P(i)
import cv2 as cv
import numpy as np
import seaborn as sns
import matplotlib.pyplot as plt
img = cv.imread('./images/kakaka.jpg')
img = cv.cvtColor(img, cv.COLOR_BGR2RGB)
img_gray = cv.cvtColor(img, cv.COLOR_RGB2GRAY)
#img_blur = cv.GaussianBlur(img_gray, (5, 5), 0)
Numpy Otsu’s thresholding
hist_cv = cv.calcHist([img_gray],[0],None,[256],[0,256]) # (256,1)
hist_np = np.bincount(img_gray.ravel(), minlength=256) # (256,)
hist_norm = hist_np/hist_np.sum()
Q = hist_norm.cumsum()
bins = np.arange(256)
fn_min = np.inf
thresh = -1
for i in range(1, 256):
p1, p2 = np.hsplit(hist_norm, [i])
q1, q2 = Q[i], 1.0 - Q[i] # Q[255] = 1.0
if q1 < 1.e-6 or q2 < 1.e-6:
continue
b1, b2 = np.hsplit(bins, [i])
# means and variances
m1, m2 = np.sum(p1*b1)/q1, np.sum(p2*b2)/q2
v1, v2 = np.sum(((b1 - m1)**2)*p1)/q1, np.sum(((b2 - m2)**2)*p2)/q2
# weighted within-class variances
fn = v1*q1 + v2*q2
if fn < fn_min:
fn_min = fn
thresh = i
print("threshold t : %d"%thresh)
threshold t : 115
mask = np.zeros_like(img_gray)
mask[np.where(img_gray >= 115)] = 1
plt.imshow(img*mask.reshape(245, 250, 1))
OpenCV otsu’s thresholding
ret, otsu = cv.threshold(img_gray, 0, 255, cv.THRESH_BINARY+cv.THRESH_OTSU)
print("threshold t : %d"%ret)
threshold t : 114
ROI
img_roi = cv.bitwise_and(img, img ,mask = otsu)