图像分割最大类间方差法

OPENCV的二值化操作中，有一种“大津阈值处理”的方法，使用函数cvThreshold(image,image2,0,255,CV_THRESH_OTSU) 实现，该函数就会使用大律法OTSU得到的全局自适应阈值来进行二值化图片，而参数中的threshold不再起作用。

OTSU算法

OTSU算法也称最大类间差法，有时也称之为大津算法，由大津于1979年提出，被认为是图像分割中阈值选取的最佳算法，计算简单，不受图像亮度和对比度的影响，因此在数字图像处理上得到了广泛的应用。它是按图像的灰度特性,将图像分成背景和前景两部分。因方差是灰度分布均匀性的一种度量,背景和前景之间的类间方差越大,说明构成图像的两部分的差别越大,当部分前景错分为背景或部分背景错分为前景都会导致两部分差别变小。因此,使类间方差最大的分割意味着错分概率最小。

设灰度图像灰度级是L，则灰度范围为[0,L-1]，利用OTSU算法计算图像的最佳阈值为：

t = Max[w0(t) * (u0(t) - u)^2 + w1(t) * (u1(t) - u)^2)]

其中的变量说明：当分割的阈值为t时，w0为背景比例，u0为背景均值，w1为前景比例，u1为前景均值，u为整幅图像的均值。

使以上表达式值最大的t，即为分割图像的最佳阈值。

以下是一段在OpenCV中实现的C语言程序，即一个使用OTSU算法提取图像阈值的函数，输入参数为一个图像指针，返回分割该图像的最佳阈值。

其中的变量说明：当分割的阈值为t时

w0为背景像素点占整幅图像的比例

u0为w0平均灰度

w1为前景像素点占整幅图像的比例

u1为w1平均灰度

u为整幅图像的平均灰度

公式:g = w0*pow((u-u0),2) + w1*pow((u-u1),2)

int MyAutoFocusDll::otsuThreshold(IplImage *frame)

{

const int GrayScale = 256;

int width = frame->width;

int height = frame->height;

int pixelCount[GrayScale];

float pixelPro[GrayScale];

int i, j, pixelSum = width * height, threshold = 0;

uchar* data = (uchar*)frame->imageData; //指向像素数据的指针

for(i = 0; i < GrayScale; i++)

{

pixelCount[i] = 0;

pixelPro[i] = 0;

}

//统计灰度级中每个像素在整幅图像中的个数

for(i = 0; i < height; i++)

for(j = 0;j < width;j++)

{

pixelCount[(int)data[i * width + j]]++; //将像素值作为计数数组的下标 }

}

//计算每个像素在整幅图像中的比例

float maxPro = 0.0;

int kk=0;

for(i = 0; i < GrayScale; i++)

{

pixelPro[i] = (float)pixelCount[i] / pixelSum;

if( pixelPro[i] > maxPro )

{

maxPro = pixelPro[i];

kk = i;

}

}

//遍历灰度级[0,255]

float w0, w1, u0tmp, u1tmp, u0, u1, u, deltaTmp, deltaMax = 0;

for(i = 0; i < GrayScale; i++) // i作为阈值

{

w0 = w1 = u0tmp = u1tmp = u0 = u1 = u = deltaTmp = 0;

for(j = 0; j < GrayScale; j++)

{

if(j <= i) //背景部分

{

w0 += pixelPro[j];

u0tmp += j * pixelPro[j];

}

else //前景部分

{

w1 += pixelPro[j];

u1tmp += j * pixelPro[j];

}

}

u0 = u0tmp / w0;

u1 = u1tmp / w1;

u = u0tmp + u1tmp;

deltaTmp = w0 * pow((u0 - u), 2) + w1 * pow((u1 - u), 2);

if(deltaTmp > deltaMax)

{

deltaMax = deltaTmp;

threshold = i;

}

}

return threshold;

}

接下来介绍OTSU方法的原理：

******************************************************************************* ******************************************************************************* *

OTSU法对于具有双峰性质的灰度图像或是彩色图像的某一通道的分割效果很好,程序为了增加健壮性加了个可以根据实际情况确定的修正值th_set.

******************************************************************************* *******************************************************************************

function y1=OTSU(image,th_set)

image=imread('color1.bmp');

gray=rgb2gray(image);%原图像的灰度图

low_high=stretchlim(gray);%增强图像，似乎也不是一定需要

gray=imadjust(gray,low_high,[]);

% subplot(224);imshow(gray);title('after adjust');

count=imhist(gray);

[r,t]=size(gray);

n=r*t;

l=256;

count=count/n;%各级灰度出现的概率

for i=2:l

if count(i)~=0

st=i-1;

break

end

end

%以上循环语句实现寻找出现概率不为0的最小灰度值

for i=l:-1:1

if count(i)~=0;

nd=i-1;

break

end

end

%实现找出出现概率不为0的最大灰度值

f=count(st+1:nd+1);

p=st;q=nd-st;%p和分别是灰度的起始和结束值

u=0;

for i=1:q;

u=u+f(i)*(p+i-1);

ua(i)=u;

end

%计算图像的平均灰度值

for i=1:q;

w(i)=sum(f(1:i));

end

%计算出选择不同k的时候，A区域的概率

d=(u*w-ua).^2./(w.*(1-w));%求出不同k值时类间方差

[y,tp]=max(d);%求出最大方差对应的灰度级

th=tp+p;

if th

th=tp+p;

else

th=th_set; %根据具体情况适当修正门限

end

y1=zeros(r,t);

for i=1:r

for j=1:t

x1(i,j)=double(gray(i,j));

end

end

for i=1:r

for j=1:t

if (x1(i,j)>th)

y1(i,j)=x1(i,j);

else

y1(i,j)=0;

end

end

end

%上面一段代码实现分割

% figure,imshow(y1);

% title('灰度门限分割的图像');

接下来介绍OTSU方法的原理以及用Ｃ语言实现：

阈值将原图像分成前景、背景两个图像。

前景：用ｎ１，ｃｓｕｍ，ｍ１来表示在当前阈值下的前景的点数，质量矩，平均灰度；背景：用ｎ２，ｓｕｍ－ｃｓｕｍ，ｍ２来表示在当前阈值下的背景的点数，质量矩，平均灰度；

当取最佳阈值时，背景应该与前景差别最大，关键在于如何选择衡量差别的标准；

而在ｏｔｓｕ算法中这个衡量差别的标准就是最大类间方差；

在以下程序中类间方差用是ｓｂ表示，最大类间方差用ｆｍａｘ，关于最大类间方差（ｏｔｓｕ）的性能：

类间方差对噪音和目标大小十分敏感，它仅对类间方差为单峰值的图像产生较好的分割效果，当目标与背景的大小比例悬殊时，类间方差准则函数可能呈现双峰或多峰，此时效果不好，但是类间方差使用是最少的。

最大类间方差（ｏｔｓｕ）的公式推导：

记ｔ为前景与背景的分割阈值，前景点数占图像比例为ｗ０，平均灰度为ｕ０；背景点数占图像比例为ｗ１，平均灰度为ｕ１；

则图像的总平均灰度为：ｕ＝ｗ０＊ｕ０＋ｗ１＊ｕ１；

前景和背景图像的方差：ｇ＝ｗ０＊（ｕ０－ｕ）．＾２＋ｗ１＊（ｕ１－ｕ）．＾２

＝ｗ０＊ｗ１＊（ｕ０－ｕ１）．＾２

此公式为方差公式。

上面的ｇ公式也就是下面程序中得ｓｂ的表达式：

当方差ｇ最大时，可以认为此时前景与背景差异最大，也就是此时的灰度值是最佳阈值；