天天看點

【OpenCV學習筆記】4.門檻值分割

  • 門檻值分割是一種基于區域的、簡單的通過灰階資訊提取形狀的技術,實作簡單、計算量小、性能穩定,是以應用廣泛。門檻值分割後的輸出圖像隻有兩種灰階值:255和0,是以門檻值分割處理又被稱為圖像二值化。

1.方法概述

  • 全局門檻值分割

    全局門檻值分割是将灰階值大于thresh(門檻值)的像素設為白色,小于或者等于thresh的像素設為黑色;或者反過來。

    假設輸入圖像為 I I I,高為 H H H,寬為 W W W, I ( r , c ) I(r,c) I(r,c)代表 I I I的第 r r r行第 c c c列的灰階值, 0 ≤ r < H , 0 ≤ c < W , 0 \leq r<H,0\leq c<W, 0≤r<H,0≤c<W,全局門檻值處理後的輸出圖像為 O O O,則:

    O ( r , c ) = { 255 , I ( r , c ) > t h r e s h 0 , I ( r , c ) ≤ t h r e s h O(r,c)= \begin{cases}255,&I(r,c)>thresh\\ 0,&I(r,c)\leq thresh \end{cases} O(r,c)={255,0,​I(r,c)>threshI(r,c)≤thresh​

    或者

    O ( r , c ) = { 0 , I ( r , c ) > t h r e s h 255 , I ( r , c ) ≤ t h r e s h O(r,c)= \begin{cases}0,&I(r,c)>thresh\\ 255,&I(r,c)\leq thresh \end{cases} O(r,c)={0,255,​I(r,c)>threshI(r,c)≤thresh​

threshold(InputArray src, OutputArray dst, double thresh, double maxval, int type)
//src-輸入矩陣,資料類型為CV_8U或者CV_32F
//dst-輸出矩陣
//thresh-門檻值
//maxval-圖像二值化時,一般為255
//type-類型
           

當type=THRESH_BINARY,

d s t ( r , c ) = { m a x v a l , s r c ( r , c ) > t h r e s h 0 , s r c ( r , c ) ≤ t h r e s h dst(r,c)= \begin{cases}maxval,&src(r,c)>thresh\\ 0,&src(r,c)\leq thresh \end{cases} dst(r,c)={maxval,0,​src(r,c)>threshsrc(r,c)≤thresh​

當type=THRESH_BINARY_INV,

d s t ( r , c ) = { 0 , s r c ( r , c ) > t h r e s h m a x v a l , s r c ( r , c ) ≤ t h r e s h dst(r,c)= \begin{cases}0,&src(r,c)>thresh\\ maxval,&src(r,c)\leq thresh \end{cases} dst(r,c)={0,maxval,​src(r,c)>threshsrc(r,c)≤thresh​

當type為THRESH_OTSU或THRESH_TRIANGLE時,src隻支援uchar類型,此時thresh作為輸出參數,即通過Otsu和TRIANGLE算法自動計算。

//examples:
threshold(src, dst, 150, 255, THRESH_BINARY);
triThe = threshold(src, dst_tri, 0, 255, THRESH_TRIANGLE+THRESH_BINARY);
           

手動設定門檻值效果不是很好,直方圖技術法、Otsu算法、熵算法是比較流行的自動選取全局門檻值的算法。

  • 局部門檻值分割

    在有些情況下,如光照不均,全局門檻值分割效果不理想,使用局部門檻值(自适應門檻值)可以産生較好的效果。

    O ( r , c ) = { 255 , I ( r , c ) > t h r e s h ( r , c ) 0 , I ( r , c ) ≤ t h r e s h ( r , c ) O(r,c)= \begin{cases}255,&I(r,c)>thresh(r,c)\\ 0,&I(r,c)\leq thresh(r,c) \end{cases} O(r,c)={255,0,​I(r,c)>thresh(r,c)I(r,c)≤thresh(r,c)​

    或者

    O ( r , c ) = { 0 , I ( r , c ) > t h r e s h ( r , c ) 255 , I ( r , c ) ≤ t h r e s h ( r , c ) O(r,c)= \begin{cases}0,&I(r,c)>thresh(r,c)\\ 255,&I(r,c)\leq thresh(r,c) \end{cases} O(r,c)={0,255,​I(r,c)>thresh(r,c)I(r,c)≤thresh(r,c)​

    局部門檻值分割的核心也是計算門檻值矩陣,常用的方法是自适應門檻值算法(移動平均值算法)。

2.直方圖技術法

一幅含有一個與背景程明顯對比的圖像具有包含雙峰的直方圖,兩個峰值對應于物體内部和外部較多數目的點,兩個峰值之間的波谷對應于物體邊緣附近相對較少數目的點。直方圖技術法就是取兩個峰值之間的波谷位置的灰階值為門檻值。下面是一種自動選取波峰和波谷的方法。

【OpenCV學習筆記】4.門檻值分割

假設輸入圖像為 I I I,高為 H H H,寬為 W W W, h i s t o g r a m I histogram_I histogramI​代表其對應的灰階直方圖, h i s t o g r a m I ( k ) histogram_I(k) histogramI​(k)代表灰階值為 k k k的像素點個數,其中 0 ≤ k ≤ 255 0\leq k\leq 255 0≤k≤255。

第一步:找到灰階值第一個峰值,即灰階直方圖最大值,并找到其對應灰階值,用 f i r s t P e a k firstPeak firstPeak表示。

第二步:找到灰階值第二個峰值,并找到其對應灰階值。第二個峰值不一定是直方圖第二大值,因為它可能出現在第一個峰值附近。可通過公式計算:

s e c o n d P e a k = a r g k m a x { ( k − f i r s t P e a k ) 2 ∗ h i s t o g r a m I ( k ) } secondPeak=arg_kmax\{(k-firstPeak)^2*histogram_I(k)\} secondPeak=argk​max{(k−firstPeak)2∗histogramI​(k)}

或 s e c o n d P e a k = a r g k m a x { ∣ k − f i r s t P e a k ∣ ∗ h i s t o g r a m I ( k ) } secondPeak=arg_kmax\{|k-firstPeak|*histogram_I(k)\} secondPeak=argk​max{∣k−firstPeak∣∗histogramI​(k)}

第三步:找到這兩個峰值之間的波谷,若出現兩個或多個波谷,取左側波谷,其對應灰階值即為門檻值。

int threshTwoPeaks(const Mat & image, Mat & thresh_out)
{
	//計算灰階直方圖
	Mat histogram = calcGrayHist(image);//函數調用
	//找到灰階直方圖最大峰值對應的灰階值
	Point firstPeakLoc;
	minMaxLoc(histogram, NULL, NULL, NULL, &firstPeakLoc);
	int firstPeak = firstPeakLoc.x;
	//尋找灰階直方圖第二個峰值對應的灰階值
	Mat measureDists = Mat::zeros(Size(256, 1), CV_32FC1);
	for (int k = 0; k < 256; k++)
	{
		int hist_k = histogram.at<int>(0, k);
		measureDists.at<float>(0, k) = pow(float(k - firstPeak), 2)*hist_k;
	}
	Point secondPeakLoc;
	minMaxLoc(measureDists, NULL, NULL, NULL, &secondPeakLoc);
	int secondPeak = secondPeakLoc.x;
	//找到兩個峰值之間最小值對應的灰階值,作為門檻值
	Point threshLoc;
	int thresh = 0;
	if (firstPeak < secondPeak)//第一個峰值在第二個峰值左側
	{
		minMaxLoc(histogram.colRange(firstPeak, secondPeak), NULL, NULL, &threshLoc);
		thresh = firstPeak + threshLoc.x + 1;
	}
	else//第一個峰值在第二個峰值右側
	{
		minMaxLoc(histogram.colRange(secondPeak, firstPeak), NULL, NULL, &threshLoc);
		thresh = secondPeak + threshLoc.x + 1;
	}
	//門檻值分割
	threshold(image, thresh_out, thresh, 255, THRESH_BINARY);
	return thresh;
}

Mat calcGrayHist(const Mat & image)
{
	Mat histogram = Mat::zeros(Size(256, 1), CV_32SC1);//存儲256個灰階級的像素個數
	int rows = image.rows;
	int cols = image.cols;//圖像的寬和高
	//計算灰階級個數
	for (int r = 0; r < rows; r++)
	{
		for (int c = 0; c < cols; c++)
		{
			int index = int(image.at<uchar>(r, c));
			histogram.at<int>(0, index) += 1;
		}
	}
	return histogram;
}
           

3.熵算法

  • 資訊熵的概念來源于資訊論,假設信源符号 u u u有 N N N種取值,記為

    u 1 , u 2 , … … , u N u_1,u_2,……,u_N u1​,u2​,……,uN​

    且每一種信源符号出現的機率,記為

    p 1 , p 2 , … … , p N p_1,p_2,……,p_N p1​,p2​,……,pN​

    那麼該信源符号的資訊熵,記為

    e n t r o p y ( u ) = − ∑ i = 1 N p i l o g p i entropy(u)=-\sum_{i=1}^{N}p_ilogp_i entropy(u)=−∑i=1N​pi​logpi​

  • 圖像也可看做一種信源,假設輸入圖像為 I I I, n o r m H i s t I normHist_I normHistI​代表歸一化的灰階圖像直方圖,那麼對于8位圖可以看作由256個灰階符号,且每個符号出現的機率為 n o r m H i s t I ( k ) normHist_I(k) normHistI​(k)組成的信源,其中 0 ≤ k ≤ 255 0\leq k\leq 255 0≤k≤255。
  • 利用熵計算門檻值步驟如下:

    第一步:計算 I I I的累加機率直方圖,又稱零階累積矩,記為

    c u m u H i s t ( k ) = ∑ i = 0 k n o r m H i s t I ( i ) , k ∈ [ 0 , 255 ] cumuHist(k)=\sum_{i=0}^{k}normHist_I(i),k\in [0,255] cumuHist(k)=∑i=0k​normHistI​(i),k∈[0,255]

    第二步:計算各個灰階級的熵,記為

    e n t r o p y ( u ) = − ∑ k = 0 t n o r m H i s t I ( k ) l o g ( n o r m H i s t I ( k ) ) , t ∈ [ 0 , 255 ] entropy(u)=-\sum_{k=0}^{t}normHist_I(k)log(normHist_I(k)),t\in [0,255] entropy(u)=−∑k=0t​normHistI​(k)log(normHistI​(k)),t∈[0,255]

    第三步:計算使 f ( t ) = f 1 ( t ) + f 2 ( t ) f(t)=f_1(t)+f_2(t) f(t)=f1​(t)+f2​(t)最大化的 t t t值,該值即為得到的門檻值,即 t h r e s h = a r g t m a x ( f ( t ) ) thresh=arg_tmax(f(t)) thresh=argt​max(f(t)),其中

    f 1 ( t ) = e n t r o p y ( t ) e n t r o p y ( 255 ) l o g ( c u m u H i s t ( t ) ) l o g ( m a x { c u m u H i s t ( 0 ) , c u m u H i s t ( 1 ) , … … , c u m u H i s t ( t ) } ) f_1(t)=\frac{entropy(t)}{entropy(255)}\frac{log(cumuHist(t))}{log(max\{cumuHist(0),cumuHist(1),……,cumuHist(t)\})} f1​(t)=entropy(255)entropy(t)​log(max{cumuHist(0),cumuHist(1),……,cumuHist(t)})log(cumuHist(t))​

    f 2 ( t ) = ( 1 − e n t r o p y ( t ) e n t r o p y ( 255 ) ) l o g ( 1 − c u m u H i s t ( t ) ) l o g ( m a x { c u m u H i s t ( t + 1 ) , … … , c u m u H i s t ( 255 ) } ) f_2(t)=(1-\frac{entropy(t)}{entropy(255)})\frac{log(1-cumuHist(t))}{log(max\{cumuHist(t+1),……,cumuHist(255)\})} f2​(t)=(1−entropy(255)entropy(t)​)log(max{cumuHist(t+1),……,cumuHist(255)})log(1−cumuHist(t))​

4.Otsu門檻值處理

  • Otsu是一種最大方差法,該算法是在判别分析最小二乘法原理的基礎上推導的。
  • 假設輸入圖像為 I I I,高為 H H H,寬為 W W W, h i s t o g r a m I histogram_I histogramI​代表其對應的灰階直方圖, h i s t o g r a m I ( k ) histogram_I(k) histogramI​(k)代表灰階值為 k k k的像素點個數,其中 0 ≤ k ≤ 255 0\leq k\leq 255 0≤k≤255。算法詳細步驟如下:

    第一步:計算灰階直方圖的零階累積矩

    z e r o C u m u M o m e n t ( k ) = ∑ i = 0 k h i s t o g r a m I ( i ) , k ∈ [ 0 , 255 ] zeroCumuMoment(k)=\sum_{i=0}^{k}histogram_I(i),k\in [0,255] zeroCumuMoment(k)=∑i=0k​histogramI​(i),k∈[0,255]

    第二步:計算灰階直方圖的一階累積矩

    o n e C u m u M o m e n t ( k ) = ∑ i = 0 k ( i ∗ h i s t o g r a m I ( i ) ) , k ∈ [ 0 , 255 ] oneCumuMoment(k)=\sum_{i=0}^{k}(i*histogram_I(i)),k\in [0,255] oneCumuMoment(k)=∑i=0k​(i∗histogramI​(i)),k∈[0,255]

    第三步:計算圖像 I I I總體的灰階平均值 m e a n mean mean,其實就是 k = 255 k=255 k=255時的一階累積矩,即

    m e a n = o n e C u m u M o m e n t ( 255 ) mean=oneCumuMoment(255) mean=oneCumuMoment(255)

    第四步:計算每一個灰階級作為門檻值時,前景區域的平均灰階、背景區域的平均灰階與整幅圖像的平均灰階的方差。對方差的衡量采用以下度量:

    σ 2 ( k ) = ( m e a n ∗ z e r o C u m u M o m e n t ( k ) − o n e C u m u M o m e n t ( k ) ) 2 z e r o C u m u M o m e n t ( k ) ∗ ( 1 − z e r o C u m u M o m e n t ( k ) ) , k ∈ [ 0 , 255 ] \sigma^2(k)=\frac{(mean*zeroCumuMoment(k)-oneCumuMoment(k))^2}{zeroCumuMoment(k)*(1-zeroCumuMoment(k))},k\in [0,255] σ2(k)=zeroCumuMoment(k)∗(1−zeroCumuMoment(k))(mean∗zeroCumuMoment(k)−oneCumuMoment(k))2​,k∈[0,255]

    第五步:找到上述最大的 σ 2 ( k ) \sigma^2(k) σ2(k),對應的 k k k即為Otsu自動選取的門檻值,即

    t h r e s h = a r g k ∈ [ 0 , 255 ] m a x ( σ 2 ( k ) ) thresh=arg_{k\in [0,255]}max(\sigma^2(k)) thresh=argk∈[0,255]​max(σ2(k))

int otsu(const Mat & image, Mat & OtsuThreshImage)
{
	//計算灰階直方圖
	Mat histogram = calcGrayHist(image);//函數調用
	//歸一化灰階直方圖
	Mat normHist;
	histogram.convertTo(normHist, CV_32FC1, 1.0 / (image.rows*image.cols), 0.0);
	//計算零階累積矩和一階累積矩
	Mat zeroCumuMoment = Mat::zeros(Size(256, 1), CV_32FC1);
	Mat oneCumuMoment = Mat::zeros(Size(256, 1), CV_32FC1);
	for (int i = 0; i < 256; i++)
	{
		if (i == 0)
		{
			zeroCumuMoment.at<float>(0, i) = normHist.at<float>(0, i);
			oneCumuMoment.at<float>(0, i) = i * normHist.at<float>(0, i);
		}
		else
		{
			zeroCumuMoment.at<float>(0, i) = zeroCumuMoment.at<float>(0, i - 1) + normHist.at<float>(0, i);
			oneCumuMoment.at<float>(0, i) = oneCumuMoment.at<float>(0, i - 1) + i * normHist.at<float>(0, i);
		}
	}
	//計算類方間差
	Mat variance = Mat::zeros(Size(256, 1), CV_32FC1);
	//總平均值
	float mean = oneCumuMoment.at<float>(0, 255);
	for (int i = 0; i < 255; i++)
	{
		if (zeroCumuMoment.at<float>(0, i) == 0 || zeroCumuMoment.at<float>(0, i) == 1)
			variance.at<float>(0, i) = 0;
		else
		{
			float cofficient = zeroCumuMoment.at<float>(0, i)*(1.0 - zeroCumuMoment.at<float>(0, i));
			variance.at<float>(0, i) = pow(mean*zeroCumuMoment.at<float>(0, i) - oneCumuMoment.at<float>(0, i), 2.0) / cofficient;
		}
	}
	//找到門檻值
	Point maxLoc;
	minMaxLoc(variance, NULL, NULL, NULL, &maxLoc);
	int otsuThresh = maxLoc.x;
	//門檻值處理
	threshold(image, OtsuThreshImage, otsuThresh, 255, THRESH_BINARY);
	return otsuThresh;
}

Mat calcGrayHist(const Mat & image)
{
	Mat histogram = Mat::zeros(Size(256, 1), CV_32SC1);//存儲256個灰階級的像素個數
	int rows = image.rows;
	int cols = image.cols;//圖像的寬和高
	//計算灰階級個數
	for (int r = 0; r < rows; r++)
	{
		for (int c = 0; c < cols; c++)
		{
			int index = int(image.at<uchar>(r, c));
			histogram.at<int>(0, index) += 1;
		}
	}
	return histogram;
}
           

5.自适應門檻值

  • 在圖像進行平滑處理時,均值平滑、高斯平滑、中值平滑用不同規則計算出以目前像素為中心的鄰域内灰階的平均值,是以可以使用平滑處理後的輸出結果作為每個像素設定門檻值的參考值。
  • 在自适應門檻值進行中,平滑算子的尺寸決定了分割出來的物體的尺寸,如果濾波器尺寸太小,那麼估計出的局部門檻值将不理想。一般平滑算子的寬度必須大于被識别物體的寬度,平滑算子尺寸越大,平滑後的結果越能更好作為每個像素門檻值的參考,當然不能無限大。
  • 假設輸入圖像為 I I I,高為 H H H,寬為 W W W,平滑算子尺寸記為 H × W H\times W H×W,其中 W W W和 H H H均為奇數,自适應門檻值分割步驟如下:

    第一步:對圖像進行平滑處理,平滑結果記為 f s m o o t h ( I ) f_{smooth}(I) fsmooth​(I),可以是均值平滑、高斯平滑、中值平滑。

    第二步:自适應門檻值矩陣 T h r e s h = ( 1 − r a t i o ) ∗ f s m o o t h ( I ) Thresh=(1-ratio)*f_{smooth}(I) Thresh=(1−ratio)∗fsmooth​(I),一般 r a t i o = 0.15 ratio=0.15 ratio=0.15。

    第三步:利用局部門檻值分割規則進行門檻值分割。

    O ( r , c ) = { 255 , I ( r , c ) > t h r e s h ( r , c ) 0 , I ( r , c ) ≤ t h r e s h ( r , c ) O(r,c)= \begin{cases}255,&I(r,c)>thresh(r,c)\\ 0,&I(r,c)\leq thresh(r,c) \end{cases} O(r,c)={255,0,​I(r,c)>thresh(r,c)I(r,c)≤thresh(r,c)​

    或者

    O ( r , c ) = { 0 , I ( r , c ) > t h r e s h ( r , c ) 255 , I ( r , c ) ≤ t h r e s h ( r , c ) O(r,c)= \begin{cases}0,&I(r,c)>thresh(r,c)\\ 255,&I(r,c)\leq thresh(r,c) \end{cases} O(r,c)={0,255,​I(r,c)>thresh(r,c)I(r,c)≤thresh(r,c)​

void adaptiveThreshold(InputArray src, OutputArray dst, double maxValue,
         int adaptiveMethod, int thresholdType, int blockSize, double C)
//src-輸入矩陣,資料類型為CV_8U
//dst-輸出矩陣
//maxValue-一般為255
//adaptiveMethod-ADAPTIVE_THRESH_MEAN_C采用均值平滑,
//               ADAPTIVE_THRESH_GAUSSIAN_C采用高斯平滑
//thresholdType-THRESH_BINARY,THRESH_BINARY_INV
//blockSize-平滑算子尺寸,奇數
//C-比例系數
           

繼續閱讀