图像分割

DSQiu

浏览: 1409900 次
性别:
来自: 广州

最近访客更多访客>>

ymgjava

spike_wq

sdsunrunner

hackxboy

博主相关

博客

微博

相册

留言

关于我

博客专栏

: NGUI所见即所得
浏览量：196932

: Effective C# ...
浏览量：0

文章分类

社区版块

存档分类

博客分类：

数字图像处理

基于边界分割基于区域分割阀值分割区域生长和分裂图像分割

终于写完数字图像分割这部分内容了，由于内容比较多，因此做一个小的内容提要，有利于更有调理的阅读，如下：

1.数字图像分割方法概要
2.基于边界分割

2.1边缘检测

2.2边界提取（简单连接，启发式搜索，曲线拟合）

3.基于区域分割

3.1阀值分割（直方图双峰，迭代法，Ostu(大律）法，基于熵的二值方法）

3.2区域生长

3.3区域分裂与合并
4.总结与实验实现（Java语言）

提供最佳阀值分割迭代算法，Otsu阀值分割算法，基于熵的二值方法的实现

最后，还是雨生的歌声，真的很不错…… 2012年9月3日2:54，更新与3：38

§1 数字图像分割处理根据不同特征，分为两类：基于边界分割和基于区域分割，主要方法有：

灰度阀值分割

边界分割法

基于纹理的分割

区域生长法

§2 基于边界分割

§2.1 边缘检测

基于边界分割其实就是点，线和边缘检测，边缘检测我在之前的一篇博文（http://dsqiu.iteye.com/blog/1638589）有介绍，可以前往查看。下面附上不同检测算子的对比：

§2.2 边界获取

使用边缘检测算子检测处理的都是孤立的一些点，需要进行边界获取将孤立点连接成线。

边缘检测的结果得到许多梯度大的点，还必须进一步处理：

原因：物体边界一般是线，而非孤立的点。必须把边缘点连接成边缘链，形成直线、曲线、各种轮廓线等，直到能表示图像中物体的边界。边界表示可以使图像的表示简洁，可以用来完成一定的文字识别任务或高层次的理解创造前提。边缘形成线特征包括二个过程：抽取可能的边缘点；将虑出的边缘连接成直线、曲线、轮廓线，或用一定的直线、曲线去拟合它们。

困难：

梯度大的点、不一定真正是边缘点。

边缘点的梯度也有可能小于周围其它点。

分岔和缺损。

对于一阶算子，需要对得到的边缘图进行细化，理想情况时，细化成单象素宽的闭合连通边界图。非理想情况下，边缘图像会有间隙存在，需要加以填充。对于二阶算子，过零点一般是单线，不需要细化，仍然需要补充间隙点。

细化算法见数学形态学。

边界跟踪

连点成线

例如，一个简单的算法：e(x,y)边缘强度，ф(x,y)边缘方向，两边缘点满足以下条件时可以连接：

|e(xi,yi)－e(xj,yj)|<T1

|ф(xi,yi)－ф(xj,yj)|mod 2π<T2

|e(xi,yi)|>T, |e(xj,yj)|>T

条件3 使小的干扰避免被误认为边缘。

从满足|e(xA,yA)|>T 的A 点出发，如果找不到满足以上条件的相邻点，算法停止。如果有多个满足条件的邻点，选边缘强度差和角度差小的点。再不断从新的起始点出发继续搜索。

启发式搜索

人工智能中的概念，是一种从多种可能的路径中选优的方法。搜索需要有评价函数，为每一条路径打分，以便于选择路径。边缘质量评价函数可以包括各点的边缘强度，也可以利用边缘的方向信息。

当有多条可能路径时，可以用全搜索（穷举法），找到评价指标最优的路径作为结果。然而，全搜索运算太大，组合爆炸。启发式搜索是在搜索的过程中设定一些启发性规则，当规则满足时就认为某一段路经合理，不再比较其它候选路经。例如：如果边缘方向上出现大的边缘点，则认为方向正确。边界的点稀疏或有长缺口时可以采用此方法。该技术对相对简单的图像效果很好，启发式搜索比全搜索减少了运算量，但不一定能找到两端点间连接的全局最优路径。

图：启发性规则的例子，只有满足图中形式的连接被认为可能。

曲线拟合

如果边缘点很稀疏，可以用分段线性或高阶样条曲线来拟合这些点，从而形成边界。

图：多边形拟合

我们都学过最小二乘法作直线拟合，更复杂的二次曲线，高斯拟合等也不难实现。数学上是某种准则下的最优参数估计问题。所用准则多为均方误差最小准则，所估计的参数是曲线方程中的未知参数。具体的拟合方法参考各种参考书。

其它边界抽取方法

Hough 变换

图搜索

动态规划

使用阀值进行图像分割

阈值分割法是一种简单的基于区域的分割技术，是一种广泛使用的图像分割技术，它利用图像中要提取的目标和背景在灰度特性上的差异，把图像视为具有不同灰度级的两类区域的组合，选取一个合适的阈值，以确定图像中每个像素点是属于目标还是属于背景。它不仅可以极大的压缩数据量，而且也大大简化了图像信息的分析和处理步骤。阈值法是首先确定一个处于图像灰度级范围内的灰度阈值T，然后将图像中每个像素的灰度值都与这个阈值T比较，根据它是否超过阈值T而将该像素归于两类中的一类。常用的方法就是设定某一阈值T，用T将图像分割成大于阈值T的像素群(目标)和小于阈值T(背景)的像素群两部分。这两类像素一般属于图像中的两类区域，所以对像素根据阈值分类达到了区域分割的目的。输入图像是F（x，y），输出图像是B（x，y），则:

§3.1 阈值化分割方法

根据图像本身的特点，可分为单阈值分割方法（全局）和多阈值分割方法（局部）;也可分为基于像素值的阈值分割方法、基于区域性质的阈值分割方法和基于坐标位置的阈值分割方法。若根据分割算法所具有的特征或准则，还可以分为直方图峰谷法、最大类空间方差法、最大墒法、模糊集法、特征空间聚类法、基于过渡区的阈值选取法等。

§3.1.1直方图阈值的双峰法

该阈值化方法的依据是图像的直方图，通过对直方图进行各种分析来实现对图像的分割。图像的直方图可以看作是像素灰度值概率分布密度函数的一个近似，设一幅图像仅包含目标和背景，那么它的直方图所代表的像素灰度值概率密度分布函数实际上就是对应目标和背景的两个单峰分布密度函数的和。图像二值化过程就是在直方图上寻找两个峰、一个谷来对一个图像进行分割，也可以通过用两级函数来近似直方图。

若灰度图像的直方图，其灰度级范围为i=0,1,…，L-1,当灰度级为k时的像素数为，则一幅图像的总像素数N为:

灰度级 i出现的概率为：

当灰度图像中画面比较简单且对象物的灰度分布比较有规律时，背景和对物象在图像的灰度值方图上各自形成一个波峰，由于每两个波峰间形成一个低谷，因而选择双峰间低谷处所对应的灰度值为阈值，可将两个区域分离。

把这种通过选取直方图阈值来分割目标和背景的方法称为直方图阈值双峰法。如下图所示，在灰度级t1和t2两处有明显的波峰，而在t处是一个谷点

具体实现的方法先做出图像f(x,y)的灰度直方图，若出现背景目标物两区域部分所对应的直方图呈双峰且有明显的谷底，则可以将谷底点所对应的灰度值作为阈值t，然后根据阈值进行分割就可以将目标从图像中分割出来。这种方法适用于适用于目标和对景的灰度差较大，直方图有明显谷底的情况。

将原始图像和阈值分割后的图像比较，可以发现有些前景图像和背景图像的灰度值太接近，导致有些前景图像没有从背景中分离出来，图像失真了。双峰法比较简单，在可能情况下常常作为首选的阈值确定方法，但是图像的灰度直方图的形状随着对象、图像输入系统，输入环境等因素的不同而千差万别，当出现波峰间的波谷平坦、各区域直方图的波形重叠等情况时，用直方图阈值难以确定阈值，必须寻求其他方法来选择适宜的阈值。

§3.1.2迭代法（最佳阀值分割迭代法）

迭代法也是聚类方法中的 K-means算法，当然下面演示的K=2的情况。

迭代式阈值选取的基本思路是:首先根据图像中物体的灰度分布情况，选取一个近似阈值作为初始阈值，一个较好的方法就是将图像的灰度均值作为初始阈值；然后通过分割图像和修改阈值的迭代过程获得认可的最佳阈值。迭代式阈值选取过程可描述如下。

(1)选取一个初始阈值T。

(2)利用阈值T把给定图像分割成两组图像，记为和。

(3)计算和均值和。

(4)选取新的阈值T,且

(5)重复第(2)~(4)步，直至和均值和不再变化为止。

具体实现时，首先根据初始开关函数将输入图逐个图像分为前景和背景，在第一遍对图像扫描结束后，平均两个积分器的值以确定一个阈值。用这个阈值控制开关再次将输入图分为前景和背景，并用做新的开关函数。如此反复迭带直到开关函数不在发生变化，此时得到的前景和背景即为最终分割结果。迭代所得的阈值分割的图像效果良好。基于迭代的阈值能区分出图像的前景和背景的主要区域所在，但在图像的细微处还没有很好的区分度。对某些特定图像，微小数据的变化却会引起分割效果的巨大改变，两者的数据只是稍有变化，但分割效果却反差极大。对于直方图双峰明显，谷底较深的图像，迭代方法可以较快地获得满意结果，但是对于直方图双峰不明显，或图像目标和背景比例差异悬殊，迭代法所选取的阈值不如其它方法。基于迭代的阈值能区分出图像的前景的主要区域所在，但在图像的细微处还没有很好的区分度。总的来说迭代法比双峰法分割的效果有很大的提高。

§3.1.3大律法(Otsu阀值分割算法）

图像记t为前景与背景的分割阈值，前景点数占图像比例为，平均灰度为；背景点数占图像比例为，平均灰度为，则图像的总平均灰度为：

从最小灰度值到最大灰度值遍历t，当t使得值

最大时的t即为分割的最佳阈值。

大津法可作如下理解：该式实际上就是类间方差值，阈值t 分割出的前景和背景两部分构成了整幅图像，而前景取值，概率为，背景取值，概率为，总均值为u，根据方差的定义即得该式。因方差是灰度分布均匀性的一种度量,方差值越大,说明构成图像的两部分差别越大,当部分目标错分为背景或部分背景错分为目标都会导致两部分差别变小，因此使类间方差最大的分割意味着错分概率最小。直接应用大津法计算量较大，因此在实现时采用了等价的公式

在测试中发现，大津法选取出来的阈值非常理想，表现较为良好。虽然它在很多情况下都不是最佳的分割，但分割质量通常都有一定的保障，可以说是最稳定的分割。

§3.1.4 类内方差最小方差法

对Otsu阀值分割做下变形，就可以得到类内最小方差法。对于每个可能的阈值，分成两类，分别计算类内方差：

使类内方差最小。

§3.1.5 最小错误概率分类法

已知前景和背景的概率分布的情况下：

例如：假设背景与前景的灰度都是正态分布，分布参数已知时可以推导出最小错误概率的门限值。假设图像上前景点的灰度值概率分布是，背景点是，选取阈值T 作为分割点。

则背景判为前景的错误：

则前景判为背景的错误：

总错误：

总错误最小，则对T 导数为零。得到：

如果前景背景都是正态分布：

是前景与背景点数量的比例，如果方差相同，比例也相同

则：

§3.1.6 基于熵的二值化方法

假设以阀值T分割图像，图像中低于阀值T的灰度的像素点构成目标物体（O），高于阀值T的像素点构成背景（B），则各个灰度级在本区的分布概率如下：

O区：

B区：
其中，L是图像灰度级总数，并且

目标和背景区域的熵分别为

对图像每个灰度级分别计算

选取使w最大的灰度级作为分割图像的阀值T。

局部自适应

当照明或透射不均匀时，整幅图像分割将没有合适的单一门限。这时，可对图像分块，对每一块分别选一门限进行分割，如果某块内有目标和背景，则直方图呈双峰。如果块内只有目标或背景，则直方图没有双峰，可根据邻居诸块得到的参数通过内插进行分割。

§3.2 基于区域的图像分割

基于边缘的图像分割：寻找区域之间的边界

基于区域的图像分割：直接创建区域

基于边缘的方法得到的结果通常不会与区域生长方法得到的分割完全一致。

二种方法结合，会是一个好办法。区域生长的方法在噪声干扰、边缘不易提取的情况下，效果更好。区域内部的一致性描述是区域生长法的基本准则。一般是其灰度，也可以考虑颜色、纹理、形状等。基于阈值的方法：基于单个点的特点。基于区域的方法考虑到相邻点的一致性。

§3.2.1 区域生长和区域合并

区域生长的原理和步骤

区域生长的基本思想是将具有相似性质的像素集合起来构成区域。具体先对每个需要分割的区域找一个种子像素作为生长的起点，然后将种子像素周围邻域中与种子像素有相同或相似性质的像素(根据某种事先确定的生长或相似准则来判定)合并到种子像素所在的区域中。将这些新像素当作新的种子像素继续进行上面的过程，直到再没有满足条件的像素可被包括进来。这样一个区域就长成了。

思路：从一些种子点生长，直到充满整个图像。种子点有监督选取。每个目标区域中至少有一个点。

需要确定：

如何选择一组能正确代表所需要区域的种子像素，如果计算结果可以看出聚类的情况，那么可以选择聚类中心作为种子像素。

生长的方法

和每次生长后的一致性准则，如灰度差小于阈值。

简单的生长方法，区域的所有8 邻域点。

如果该点加入后满足一致性准则，则加入。两个区域满足一定准则，可以合并，该准则可以考虑两个区域分别的均值和方差。如果没有预先确定的种子点，可采用一般步骤

1. 用某种方法把图像分割成许多小区域。

2. 定义合并相邻区域的准则。

3. 按照合并准则合并所有相邻的区域，如果没有再能够合并的块后停止。

采用不同的初始分割方法和合并准则，可以得到适应不同情况的算法。另外，区域合并得结果通常还依赖于区域合并的顺序。

相邻区域的特征值之间的差异是计算边界强度的一个尺度。如果给定边界两侧的特征值差异明显，那么这个边界很强，反之则弱。强边界允许继续存在，而弱边界被消除，相邻区域被合并。计算是一个迭代过程，每一步重新计算被扩大的区

域成员隶属关系，并消除弱边界。没有弱边界可消除时，合并过程结束。

过程看起来象一个物体内部区域不断增长，直到到达边界为止的过程。

该方法计算开销大，但能够同时利用图像的若干种性质（多种描述），对自然景物分割方面效果相对最优。

生长准则和过程

区域生长的一个关键是选择合适的生长或相似准则，大部分区域生长准则使用图像的局部性质。生长准则可根据不同原则制定，而使用不同的生长准则会影响区域生长的过程。下面介绍3种基本的生长准则和方法。

(1) 灰度差准则

区域生长方法将图像的像素为基本单位来进行操作，基本的区域灰度差方法主要有如下步骤:

①对图像进行逐步扫描，找出尚没有归属的像素；

②以该像素为中心检查它的邻域像素，将邻域中的像素逐个与它比较，如果灰度差小于预先确定的值，将它们合并；

③以新合并的像素为中心，返回到步骤②，检查新像素的邻域，直到区域能进一步扩张；

④返回到步骤①，继续扫描直到不能发现没有归属的像素，则结束整个生长过程。

采用上述方法得到的结果对区域生长起点的选择有较大依赖性。为克服这个问题可采用下面的改进方法:

①设灰度差的阈值为零，用上述方法进行区域扩张，使灰度相同像素合并；

②求出所有邻接区域之间的平均灰度差，并合并具有最小灰度差的邻接区域；

③设定终止准则，通过反复进行上述步骤②中的操作将区域依次合并直到终止准则满足为止。

另外，当图像中存在缓慢变化的区域时，上述方法有能会将不同区域逐步合并而产生错误。为克服这个问题，可不用新像素的灰度值去与邻域像素灰度值比较，而用新像素所在区域的平均灰度值去与各邻域像素的灰度值进行比较。

对一个含N个像素的图像区域R，其均值为:

对像素是否合并的比较测试表示为:

其中T为给定的阈值。

区域生长的过程中，要求图像的同一区域的灰度值变化尽可能小，而不同的区域之间，灰度差尽可能大。下面分两种情况进行讨论:

1)设区域为均匀的，各像素灰度值为m与一个零均值高斯噪音的叠加。当测试某个像素是否合并时，条件不成立的概率为:

这就是误差概率函数，当T取3倍的方差时，误判概率为1-99.7%。这表明，当考虑灰度均值时，区域内的灰度变化应尽量小。

2)设区域为非均匀，且由两部分不同目标的图像像素构成。这两部分像素在R中所占比例分别为 ,灰度值分别为，则区域均值为。对灰度值为m的像素，它与区域均值的差为：

可知正确的判决概率为:

这表明，当考虑灰度均值时，不同部分像素间的灰度差距应尽量大。

(2) 灰度分布统计准则

这里考虑以灰度分布相似性作为生长准则来决定区域的合并，具体步骤为:

①把图像分成互不重叠的小区域；

②比较邻接区域的累积灰度直方图根据灰度分布的相似性进行区域合并；

③设定终止准则，通过反复进行步骤②中的操作将各个区域依次合并直到终止准则满足。

这里对灰度分布的相似性常用两种方法检测(设分别为两邻接区域的累积灰度直方图):

Kolmogorov-Smirnov检测:

Smoothed-Difference检测:

如果检测结果小于给定的阈值，即将两区域合并。

采用灰度分布相似合并法生成区域的效果与微区域的大小和阈值的选取关系密切，一般说来微区域太大，会造成因过渡合并而漏分区域:反之，则因合并不足而割断区域。而且，图像的复杂程度，原图像生成状况的不同，对上述参数的选择有很大影响。通常，微区域大小q和阈值T由特定条件下的区域生成效果确定。

(3) 区域形状准则

在决定对区域的合并时也可以利用对目标形状的检测结果，常用的方法有两种:

①把图像分割成灰度固定的区域，设两邻区域的周长分别为p1和p2 ，把两区域共同边界线两侧灰度差小于给定值的那部分长度设为L，如果(t1,为预定阈值):

则合并两区域;

②把图像分割成灰度固定的区域，设两邻接区域的共同边界长度为B，把两区域共同边界线两侧灰度差小于给定值的那部分长度设为L，如果(为预定阈值）

则合并两区域。

上述两种方法的区别是:第一种方法是合并两邻区域的共同边界中对比度较低部分占整个区域边界份额较大的区域，而第二种方法则是合并两邻接区域的共同边界中对比度较低部分比较多的区域。

§3.2.2 区域分裂

区域分裂是与区域合并相反的过程。

先假设整个图像是一个对象。

不满足一致性准则，则分裂，（一般是分裂成4 个子图像）。

分裂过程反复进行，直到所有分开的区域满足一致性准则。

如果图像是大小为N×N 的正方形，N 是2 的整数次幂，即N＝2n，则分裂得到的所有区域都是2 的整数次幂的正方形。逐层得到的一种图像的表示方法，叫做四叉树（quadtree）。每个节点有4 个子节点。

四叉树是一种很方便的区域描述方法。

分裂的结果，可能有相邻且满足一致性准则，但分到了不同的块中。解决方法―增加合并过程。

进一步：分裂/合并综合的方法。

§3.2.3 分裂和合并法

如果用树表示一幅图像，树根代表整个图像，树叶代表每个象素。

分裂的方法是从树根开始处理；合并的方法是从树叶开始处理。

如果中间层次开始处理，按照一致性准则，该合并的合并，该分裂的分裂，无疑是更好的选择。

减少计算量的同时，具有分裂和合并法的优点。

方法的要素：输入图像、一致性准则、初始的区域

划分

以四叉树的某一层节点作初始区域划分。

1. 合并具有一致属性的共根四个节点块。

2. 分开不满足一致性要求的块，上一步没有合并得块，如果它的四个子块不满足一致性，则将其分解成4 个子块。分出的子块如不满足一致性，还可以继续分解。

3. 对相邻的具有一致属性的块合并，即使不在同一层或者没有共同的父节点。将四叉树变为邻接图表示。

4. 如果没有进一步的合并或分解，算法终止。

§3.4 边缘与区域相结合的分割

边缘与区域组合分割的主要思想是结合二者的优点，通过边缘点的限制，避免区域的过分割；同时，通过区域分割补充漏检的边缘，使轮廓更加完整。例如：先进行边缘检测与连接，再比较相邻区域的特征(灰度均值、方差等)，若相近则合并；对原始图像分别进行边缘检测和区域生长，获得边缘图和区域片段图后，再按一定的准则融合，得到最终分割结果。

连通区域标记

图像分割一般得到的多个区域，通常需要通过标记每个象素分别属于哪个区域，分别把每个区域提取出来。

在有多个对象需要处理的情况下，这一步骤是必不可少的。

基于边缘的方法已经得到闭合边界，可以采用边界跟踪和内部填充的方法。见platt 书P581。

基于区域的方法一般采用连通性分析方法，按照一定的顺序把连通的像素用相同的序号标注。

算法：

1. 把所有像素点放到待处理点集合A 中

2. 如果 A 空则结束。否则从A 中任意移出一点,作为连通域a（用集合表示）的初始点。

3. 在 A 中寻找所有与a 连通的点，并移到a 中。如果没有找到，转到2，寻找下一个连通域

4. 转到 3，迭代寻找新的连同点。算法逻辑简单，速度慢。快速算法：每点只需要遍历一次。

图像分割的评价

找到某种方式来评价图像分割方法的好坏对选择和研究算法有非常有用。

Haralick 和Shapiro 建立了下列定性的指导：

对于某些特征（如灰度、纹理），同一个区域的图像应该一致和均匀

区域内部应该简单，没有很多空洞

相邻的区域在满足区域内部一致性的特征上应该有显著的区别

每个区域的边界应该简单而不粗糙，并且空间位置准确

目前仍没有定量的度量标准。

§4.1图像分割总结

算法的其他的分类方式：

本文以使用特征的物理意义不同来区分，有利于掌握本质特点。

章毓晋书按照计算过程分：串行、并行。并行算法所有判断可以独立同时做出，串行算法后一步要利用前一步的结果。

总结

1. 图像分割是图像描述的重要基础。

2. 怎样获得分割呢？

可以设想分割后的区域对某种或某些图像特性具有可靠的一致或平稳；相邻区域之间的边界应该完整、不破碎又是能精确定位的。

3. 基于边界或区域的方法都是基于基本假设，假设不成立时就会遇到困难。基于边界的方法，物体内部边缘，光照等原因物体间边缘差别小或平滑过渡。人可以在理解的基础上判别，需要更上层知识的指导。基于区域的方法，区域内部特征不一定那么均匀。

4. 虽然图像处理一开始就研究分割的问题，也取得了相当进展，但尚无一种适合于所有图像的通用分割算法。

5. 人们至今仍然在努力发展新的、更有潜力的方法，期待更通用、更完美的分割结果。

6. 本章介绍了一些典型的基本的分割算法，实际上，现有的方法通常是针对具体问题。为了得到好的性能，要根据实际问题选择或设计算法。众多方法的关键是模型和条件，不同的模型和条件决定了不同的方法。根据具体任务找到或设计匹配的模型，根据条件求解。

7. 在研究新方法的同时，分割参数的自动选择，利用统计理论作决策，适用于快速处理的数据结构，较为通用的分割结果质量评价标准都是研究者关心的内容。如果区域内部灰度有剧烈的变化，则需要用纹理分析来分割。

§4.2 程序实现

一维最大熵分割算法

//一维最大熵分割算法
    public int segment(int[] pix, int w, int h)
	{
		int i, j, t;
        double a1, a2, max, pt;
        double[] p = new double[256];
        double[] num = new double[256];	    
	       
  	    int[][] im = new int[w][h];
		
		for(j = 0; j < h; j++)
		    for(i = 0; i < w; i++)		    
				im[i][j] = pix[i+j*w]&0xff;
  	    
        for (i = 0; i < 256; i++)
            p[i] = 0;
        
        //统计各灰度级出现的次数
        for (j = 0; j < h; j++)
            for (i = 0; i < w; i++)            
                p[im[i][j]]++;
        
        /**
         *关于计算各灰度级出现的概率
         *1.因为(p[j]/(w*h)) / (pt/(w*h)) = p[j] / pt
         *  所以计算p[j] / pt不必计算概率
         *2.当p[j]=0时，计算Math.log(p[j] / pt)将出现无穷大.但
         *  此时p[j] / pt) * Math.log(p[j] / pt)=0
         *  所以在计算a1时,不必计算这一项       
         */
        int hw =  h*w;               
        for (i = 0; i < 256; i++)
        {
            a1 = a2 = 0.0;
            pt = 0.0;
            for (j = 0; j <= i; j++)
                pt += p[j];
            
            for (j = 0; j <= i; j++)
            
            	if(p[j]>0)
                    a1 += (p[j]/pt) * Math.log(pt/p[j]);
                        
            for (j = i+1; j <256; j++)
              	if(p[j]>0)
                    a2 += (p[j] /(hw-pt))* Math.log((hw - pt)/p[j]);
            
            num[i] = a1 + a2;
        }
        
        max = 0.0; t = 0;
        for (i = 0; i < 256; i++)
        {
            if (max < num[i])
            {
                max = num[i];
                t = i;
            }
        }        
        return t;
	}

二维最大熵分割算法, 使用递推算法

   public int segment2(int[] pix, int w, int h)
	{
		int i, j, u, v, t;
        double a1, a2, max, pa, pb, pa2, pb2, sum;
        double[][] p = new double[256][256];
        double[][] num = new double[256][256];
	    
	    int[][] im = new int[w][h];
		
		for(j = 0; j < h; j++)
		    for(i = 0; i < w; i++)		    
				im[i][j] = pix[i+j*w]&0xff;
				
        for(i = 0; i < 256; i++)
            for(j = 0; j < 256; j++)
                p[i][j] = 0;
        
        //统计2维直方图p[i][j]
        for(j = 1; j < h-1; j++)
        {
            for(i = 1; i < w-1; i++)
            {
            	t = (int)((im[i-1][j]+im[i+1][j]+im[i][j-1]
            	  +im[i][j+1]+im[i][j])/5);//4-邻域均值
                p[im[i][j]][t]++;
            }
        }
        
        pa = 0.0; pb = 0.0; max = 0.0; t = 0;
        for(i = 49; i < 200; i=i+2)
        {   
            System.out.println((int)(i*100/199)+" %");
            for(j = 0; j < 256; j++)
        	{
        		a1 = 0.0; a2 = 0.0; 	            
	            pb = 0.0;
	            
	            //递推算法计算pa
	            if(j != 0)
	            {
	            	for(u = 0; u <= i; u++) 
	                    pa += p[u][j];
	            }
	            else
	            {   
	                pa = 0.0;
	            	for( u = 0; u <= i; u++)
	            	    pa += p[u][0];
	            }
	            	            
	            //递推算法计算pb
	            if(j != 0)        
	            {
	            	for(u = i+1;u < 256;u++)
	            	    pb -= p[u][j];
	            }
	            else
	            {
	            	pb = 0;
	            	for(u = i+1;u < 256;u++)
	            	    for(v = j+1; v < 256; v++)
	            	        pb += p[u][v];
	            }
	           
	            for(u = 0; u <= i; u++)
	               	for(v = 0; v <= j; v++)
	               	    if(p[u][v] > 0)
	                        a1 += (double)(-p[u][v]/pa)* Math.log(p[u][v]/pa);
	           
	            for(u = i+1; u < 256; u++)
	               	for(v = j+1; v < 256; v++)
	               	    if(p[u][v] > 0)
	               	        a2 += (double)(-p[u][v]/pb)* Math.log(p[u][v]/pb);  
	           
	            num[i][j] = a1 + a2;	                        
            }
        }
        
        max = 0.0; t = 0;
        for (i = 0; i < 256; i++)
        {
        	for(j = 0; j < 256; j++)
        	{
                if (max < num[i][j])
                {
                    max = num[i][j];
                    t = i; 
                }
            }
        }    
        return t;
	}

最佳阈值分割

	public int bestThresh(int[] pix, int w, int h)
	{
		int i, j, t,
		    thresh, 
		    newthresh,
		    gmax, gmin;         //最大,最小灰度值
        double a1, a2, max, pt;
        double[] p = new double[256];
        long[] num = new long[256];
	
	    int[][] im = new int[w][h];
		
		for(j = 0; j < h; j++)
		    for(i = 0; i < w; i++)		    
				im[i][j] = pix[i+j*w]&0xff;
				
        for (i = 0; i < 256; i++)
            p[i] = 0;
        
        //1.统计各灰度级出现的次数、灰度最大和最小值
        gmax = 0;
        gmin =255;
        for (j = 0; j < h; j++)
        {
            for (i = 0; i < w; i++)
            {
            	int g = im[i][j];
                p[g]++;
                if(g > gmax) gmax = g;
                if(g < gmin) gmin = g;
            }
        }
        
        thresh = 0;
        newthresh = (gmax+gmin)/2;
        
        int meangray1,meangray2;
        long p1, p2, s1, s2;
        for(i = 0; (thresh!=newthresh)&&(i<100);i++)
        {
        	thresh = newthresh;
        	p1 = 0; p2 = 0; s1 = 0; s2 = 0;
        	
        	//2. 求两个区域的灰度平均值
        	for(j = gmin; j < thresh;j++)
        	{
        		p1 += p[j]*j;
        		s1 += p[j];        		
        	}
        	meangray1 = (int)(p1/s1);
        	
        	for(j = thresh+1; j < gmax; j++)
        	{
        		p2 += p[j]*j;
        		s2 += p[j];        		
        	}
        	meangray2 = (int)(p2/s2);
        	//3. 计算新阈值
        	newthresh = (meangray1+meangray2)/2; 	
        }
        return newthresh;
	}

Otsu阈值分割

	public int otsuThresh(int[] pix, int iw, int ih)
	{
		ColorModel cm = ColorModel.getRGBdefault();
        int wh = iw * ih;
        int[][] inIm = new int[iw][ih]; 
 
        int i, j, t;
        int L = 256;
        double[] p = new double[L];
                        
        for (j = 0; j < ih; j++)
            for (i = 0; i < iw; i++)
                inIm[i][j] = pix[i+j*iw]&0xff;               

        for (i = 0; i < L; i++)
            p[i] = 0;

        //计算各灰度出现次数
        for (j = 0; j < ih; j++)
            for (i = 0; i < iw; i++)
                p[inIm[i][j]]++;

        //计算各灰度级出现概率
        for (int m = 0; m < L; m++)
            p[m] = p[m] / wh;

        double[] sigma = new double[L];
        for (t = 0; t < L; t++)
        {
            double w0 = 0;
            for (int m = 0; m < t+1; m++)
                w0 += p[m];
            double w1 = 1 - w0;

            double u0 = 0;
            for (int m = 0; m < t + 1; m++)
                u0 += m * p[m] / w0;

            double u1 = 0;
            for (int m = t; m < L; m++)
                u1 += m * p[m] / w1;

            sigma[t] = w0*w1*(u0-u1)*(u0-u1);
        }
        double max = 0.0;
        int T = 0;
        for (i = 0; i < L-1; i++)
        {
            if (max < sigma[i])
            {
                max = sigma[i];
                T = i;
            }
        }        
        return T;        		
	}

上述算法函数返回值都是图像分割的阀值，然后进一步对图像进行阀值灰度处理。

 //图像序列pix阈值分割	
	public int[] thSegment(int[] pix, int iw, int ih, int th)
	{						
		int[] im = new int[iw*ih];
		int t;
		for(int i = 0; i < iw*ih; i++)	
		{
			t = pix[i]&0xff;
								
			if(t > th) 
				im[i] = (255<<24)|(255<<16)|(255<<8)|255;//背景色
			else
			    im[i] = (255<<24)|(0<<16)|(0<<8)|0;      //前景色为		
		}
		return im;
	}

最后附上几张实验效果图：