阈值确定方法

图形（或图像）在计算机里主要有两种存储和表示方法。矢量图是使用点、 直线或

多边形等基于数学方程的几何对象来描述图形， 位图则使用像素来描述图 像。一般来说，照片等相对杂乱的图像使用位图格式较为合适， 矢量图则多用于 工程制图、标志、字体等场合。矢量图可以任意放缩，图形不会有任何改变。而 位图一旦放大后会产生较为明显的模糊，线条也会出现锯齿边缘等现象。

矢量图从本质上只是使用曲线方程对图形进行的精确描述，

在以像素为基本

显示单元的显示器或打印机上是无法直接表现的。 将矢量图转换成以像素点阵来 表示的信息，再加以显示或打印，这个过程称之为栅格化（ Rasterization），见图 1。

栅格化的逆过程相对比较困难。假设有一个形状较为简单的图标，保存成一 定分辨率的位图文件。我们希望将其矢量化，请你建立合理的数学模型，尽量准 确地提取出图案的边界线条，并将其用方程表示出来。

二、 问题分析

本题的要求是完成位图的矢量化，通过建立合理的数学模型，将一个有一定 分辨率

的位图文件尽量准确地提取出图案的边界线条， 最终将位图用方程的形式 表示出来。解决本问题的流程图见下 图。首先，通过MATLAB读取位图的各个 像素的像素值（0-1）,得到位图各个点的灰度值，通过最大类间方差法和最大熵 法确定阈值，完成灰度的二值化，使各个像素点的灰度值全部由0或1表示。其 次，将位图的轮廓通过合适的算法提取出来，根据特征值对轮廓进行拟合。最后， 根据拟合的函数完成位图的矢量图， 完成其矢量化过程，并通过对比矢量图和原 始位图对应的。

二、问题假设及符号说明

3.1问题假设 3.2符号说明

四、模型建立

4.1模型准备

本题要求将一个形状较为简单的图标， 保存成一定分辨率的位图文件，即将 位图矢量化。

阈值：指释放一个行为反应所需要的最小刺激强度， 值化的临界值。 4.2阈值的确定方法 421最大类间方差法

本文指像素点灰度值二

最大类间方差法的基本思想是将待分割图像看作是由两类组成的整体， 一类 是背景，一类是目标 ⑹。用各个像素点像素值的方差来衡量目标和背景之间额 差别。方差越大，即两像素值的区别越明显，不能分在同一类别中。本问中，为 使得目标和背景两类的类间方差最大的灰度级别即认为是最佳阈值。 一副大小为 M N个像素、灰度级别为L的图像，若设图像中灰度级为i的像素个数为Ni， 则灰度级i的概率为

P

N

i

N M

根据最大类间方差法的分割思想，最佳阈值公式为

* 2 2

T arg max 巳

0 t L 1

w° Pb wb w°

其中，

T

Pa P

i 1

L

Pb P

i T 1

Pa

L

Wb

P i —L

i T 1 P

b

L

Wo iPi

i 1

一维最大类间方差法是建立在一维灰度值分布图基础上的确定阈值的方法， 而一维灰度值分布图没有考虑图像的空间信息， 仅仅表示出了各个像素点的 灰度分布，但是图像中的每个像素与其相邻的像素是有相关性的，

所以在以为灰

度值分布图基础之上的阈值分割方法所建立的分割模型考虑的因素不全面。因 此，在对图像进行分割时，容易受到噪声的干扰，得到的分割结果和阈值结果不 精确，将影响后

续位图矢量化的准确度。

为了增强算法的抗噪性能、提高其确定阈值结果的精确度，本文将最大类间 方差法从一维算法推广至二维算法。二维最大类间方差法是建立在二维灰度值分 布图的基础上计算阈值的算法，下面对其进行介绍。

设一幅大小为M N个像素的图像灰度级为L， f x, y为图像在x, y的灰 度值，g x, y是以x, y为中心，k k邻域内的平均灰度值，其表达式为

k 2

g x, y

k

2

m,y

k 2

式中，k为邻域的大小;

g x, y为邻域内的平均灰度值;

f x m, y n为图像在点x m, y n的灰度值;

M为图像的宽度； N为图像的高度。

f x, y和g x, y即形成了一个二元组f x, y , g x, y ，令Nij为图像中 f x,y i且g x, y j的像素的个数0 i,j L 1 ,则其联合概率密度为

Nj

Pj

M N

j的像素联合概率密度;

式中，Rj为图像中f x, y i且g x, y

Nij为图像中f x, y i且g x, y j的像素的个数。

由最大类间方差法求出的阈值向量 s,t会将直方图分为四个部分，分别为目 标部分、背景部分、边界部分、噪声部分。在大多数情况下，近似的认为边界部 分和噪声部分的概率为0,因此这种近似在多数情况下是合理的。

图像目标出现的概率为

Ro

图像背景出现的概率为

L L

R

R

i s 1j t 1

式中，Rj为图像中f x, y i且g x, y

L为图像灰度级

j的像素联合概率密度;

图像总均值矢量为

T TO， T1

iPj, jRij

图像目标和背景额均值矢量为

st

O

OO ， 01

i 1 j 1

st

iRj,

i 1 j 1

jPij

1 10， 11

iPij,

i s 1j t 1 i s 1j t 1

jPij

因为假设边界区域和噪声区域的概率为 0,因此有

Po P1 1 P

T o 0

p

1 1

式中，o和1分别为目标区域和背景区域的均值矢量

类间方差定义为

SB

p

o

p

1 1 T 1 T

取类间方差的迹作为测度，有

2 2

p

trSB

式中，

s t

i 0 T0 j

p

i T1

P0 1 P0

jRj

i s 1j t 1

式中，Rj为图像中f x,y i且g x, y j的像素联合概率密度;

L为图像灰度级。 最佳阈值公式为 s*,t* argmaxtrSB s,t

1 s L 1 t L

式中，trSB为测度。

遍历图像的灰度级，求使目标和背景的类间方差的迹最大的灰度级， 佳分割阈值s*,t* 。 4.2.2最大熵法

由熵的定义可知，当一个信息源中所有的事件以相同的概率出现时， 这

时候信息熵大，因为信息源的不确定性大，所以信息熵大。用最大熵法求最佳阈 值，就是求一个分割阈值使得目标和背景两类的信息熵之和最大，一幅大小为

M N个像素、灰度级为L的图像，设图像中灰度级为i的像素个数为Ni，则灰 度级i的概率为

Ni N M

图像目标区域的熵和图像背景区域的熵可以表示为

即得最

Hf

Hb

log 1 Pt

P

式中，Hf为图像目标区域的熵;

H b为图像北京区域的熵;

Pt

最大熵法的最佳阈值公式为

*

T arg max H f T Hb T

遍历图像的灰度级，求目标区域的熵和背景区域的熵的和最大的灰度级，即 得最佳分割阈值T* o

424确定阈值方法总结

将以上提出的累积剩余熵阈值分割算法用于实际图像分割， 将分割结果与经 典的最大类间方差法和最大熵法进行了比较。 实验中用了人、直升机、车、帆船、 船和发电站图。表1是对分割结果的主观评价，表2是三种阈值分割方法的最佳 分割阈值的比较。

表1分割结果主观评价表

图像 人 直升机 车 帆船 船 发电站

最大类间方差法

差 差 差 一般 差 差

最大熵法

好 好 好 一般 一般 一般

最大累积剩余熵法

好 好 好 一般 好 好

表2最佳分割阈值比较表

图像 人 直升机 车 帆船 船

最大类间方差法

129 82 65 134 99

最大熵法 165 94 85 156 126

最大累积剩余熵法

153 94 77 137 144

由上表分析可知，图像的不同对应的最佳分割阈值的确定方法不同， 女口：

当

图像为人、直升机或者汽车时，最大熵和最大累积剩余熵确定阈值的方法比最大 类间方差法结果更精确。 4.3二值化

根据像素点的灰度值f x,y和阈值T ,可将图像二值化。

fx,y

1

“ 0 f x,y T

T

式中，f x, y为像素的灰度值，T为最佳分割阈值 4.4基于细化算法的轮廓提取 4.4.1基础术语

为简便的描述图像中像素之间的相互关系，本文引用了邻域和连通的概念 [8]。

・邻域和4相邻

对于图像中的某个像素点p m,n，则与之在水平方向（左和右）和垂直方

向（上和下）相邻的4个像素点坐标分别为 m 1,n，m 1, n，m, n 1，m,n 1 ， 则这4个像素点组成了像素点p m,n的4邻域，表示为N4 p。而这4个像素点 在位置上就与像素点p m,n相邻

• (m,rU 1) • O • (m-1 ,n) (m?n) (m+lji) •

图4邻示意图 图4邻坐标关系图

•8邻域和8相邻

若取像素p m,n四周的8个像素点作为相邻点，则 像素点p m, n的这8个

(ni-Lii-l) (m.n-1) 相邻点就构成了 8邻域，用N8 p表示。

• • • • • O • • • 图8邻坐标关系图

图8邻示意图

•像素间的邻接

像素的相邻仅说明了两个像素在位置上的关系，若再加上取值相同或相近, 则称两个像素邻接。

① 两个像素点p和q邻接的条件

像素点p m,n和q s,t位置上满足相邻，即 4 相邻：m, n N4 q 或者 s,t N4 p 8 相邻：m, n N8 q 或者 s,t N8 p

且像素点p m,n和q s,t灰度值相近，即称为灰度值相近（似）准则，即

p V和q V，其中V

② 邻接的类别 4邻接：p, q V且s,t

4

Vi\",…

N4 p，则称像素点 p m, n和q s, t 4邻接，记为

p q ；

8邻接：p,q V且s, t p，则称像素点 p m, n和q s,t 8邻接，记为

8

p q。

两种邻接的关系：4邻接必8邻接，反之不一定成立。

两种邻接及其关系见下图，相似性准则为 V 1 ,其中像素点p和q4邻接，也

8邻接，点q和r 8邻接但非4邻接

P 0 0 1 1 0 0 0 0 0 0 1 1 0 0 0 0 / 1| 0 0 0 q 1 0 0 ,1 0

r ⑴像素标记

1 （2）像素取值

图邻接关系图

⑶4邻接 ⑷8邻接

•像素间的连通 ① 通路

设p m,n与q s,t之间的各像素点形成的连线L为

L p.q m,n m°,n° , m^ m ，…,mij 1, ni 1 , mi, nj，…，m.,n. s,t

若 mj i,nj i与 m,nj邻接1 i n，贝U L p,q 称为p m,n与q s,t之间的一 条通路，N为通路长度。与链接一样，通路也分为 ② 连通性

若S是图像中的一个子集，p,q S，且存在一条由S中像素组成的从p到q 的通路，则称p在图像集S中与q连通，连通也分为4连通和8连通。

4通路和8通路。

0 0 0 0 0 0 0 0 0 0 1 1 0 d 0 ■ 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 1T 1 0 1 0 0 0 0 1 1 0 0 0 0 1 0 1\\ 0 0 S- 0 0 1- -1- -1 0 0 0 0 0 1 0 1 1 0 0 0 1 0 八/ - 0 -s 0 0 0 0 0 1 0 1- -1 -1 0 0 - 0 0 0 0 0 0 0 0 0 0 0 0 0 0 （坊8连通

⑻4连通 其中尸⑴

图连通示意图

442二值图像的细化

所谓细化，就是在保持原图像拓扑结构的情况下， 将线性图转成单像素宽的 骨架。细化算法是一种特殊的多次收缩的过程， 由于它在收缩过程中不会消去只 有一个临接点的边界点，因此可以保证收缩过程中的连通性，细化后的曲线位于 区域，一般采取由左至右、由上至下的细化过程，反复多次地消去具有两个 以上临接点的边界点。其基本思想是对某一像素p，观察其邻点像素的分布情况， 在不破坏图形连通性的情况下，将 p点的值由1变为0。下图为细化过程实例。

i I 1 I I I

I

i

I

1 i i (a)原始图像 (b) 消去上边界后的区域图像

(c)

消去下边界后的区域图像 (d) 消去左边界后的区域图像 图细化过程示例图

(e)消去右边界后的区域图像

4.5拟合曲线

通过对图像的细化得到一个单像素宽的骨架轮廓， 系列然而，若将组成骨架的一 点全部拟合到一个函数中，拟合结果的精确次数将会过高， 不便于计算机的 运算。因此，为简化运算、降低函数的次幕，我们对构成单像

素宽骨架的所有点 进行了基于特征的分类和提取，为保证结果的简洁性，根据不同类型点的特征对

骨架进行分割，分别拟合各段轮廓，得到幕级数较低的拟合函数曲线 4.5.1提取分类点

根据构成单像素宽骨架的点的特征，将点分为三类：间断点、分割点、特征 点。 •可断点

间断点为轮廓中的突然出现的尖点， 类似于数学概念中的不可导点，由于图 像中的像素点都是离散的，因此对此类点也进行离散化处理。本文引用差分的方 法，用一阶差商来近似一阶偏导数，一阶差商有三种形式：一阶前差商、一阶后 差、一阶中心差商，分别记为

p Pi i Pi x x P Pi Pi i x x

■分割点

分割点为轮廓中曲率突变的点，类似于数学概念中刻画离散点的曲率， 可用 于将骨

架轮廓分割成段。曲线中变化剧烈的点的曲率较大，在保证精度的前提下， 变化曲率较大的地方要特殊保留，可以通过估算曲线上各点的余弦函数值来间接 地估算该点的曲率，以此得到曲线上变化剧烈的点作为明显分割点进行先行保 留。

设P点是轮廓曲线连通段上点序中第i个离散点，Pi k和Pi k分别是距离Pi 点步长为k的前后两个点，贝U Pi点的曲率值可用下列公式估算。

cos Xi / 2 2 / 2

：X 人 k X k 人，％

x k Xi k x % w k w k % % k Wk

2

yi

根据阈值大小，将曲率值 大于阈值的点作为分割点。 •特征点

特征点为轮廓中能最大限度代表整个轮廓形状的点，类似于近似轮廓组成 点。为避免拟合函数中出现有高度波动性的高次多项式， 特征点的选取可简化轮 廓，最终降低拟合函数的次数。

以图为例进行拟合，特征点的选取要经历如下步骤。

确定一个拟合误差 ，连接始点A和终点H，找出距离AH直线的最远点C， 若C到AH的距离小于，则AH即为满足要求的拟合矢量，若C到AH的距离 大于，则C到AH之间的一个锚点，可用 AC和CH两个矢量代替AH矢量。 用同样的方法依次分别对 AC和CH处理，直到所有离散点与其拟合矢量的距离 均小于预定的误差为止，最终以上满足预定误差 条件的离散点即为单像素宽

轮廓的特征点

D

图 基于距离判断的特征点提取示例

用有限次端点迭代的方法，找出单像素宽轮廓图中满足一定精度的拟合矢 量。以图为例，具体的取特征点并矢量化的过程如下

[7] :（ 1

16）

①以AH为拟合矢量，找出距离 AH的最远点C，C距AH> ，贝U C为锚点; ② 以AC、CH代替AH，先找出AC区段内所有点位距AC的最远点位B，B距 AC> ，则B仍为锚点；

③ 以AB、BC代替AC ,先找出AB区段内距AB最远点位并发现该点位距 AH> 。 则B点为一个顶点。再找BC区段内距BC的最远点位并发现该点距 BC< ，则 C点为顶点。至此AC区段已处理完毕；

④ 找出CH区段内距CH的最远点位D，D距CH> ，贝U D为锚点；

⑤ 以CD、DH歹徒CH，先找CD区段内距CD的最远点位并发现该点距 CD< ， 则D为顶点。再找DH区段内距DH最远点位E，E距DH> ，贝U E为锚点；

⑥ 以DE、EH代替DH，先找DE区段内距DE的最远距离点位并发现该点距 DE< ，贝U E为顶点。再找EH区段内距EH的最远距离点F, F距EH> ，F为 锚点；

⑦ 以EF、FH代替EH，先找EF区段内距EF的最远点位并发现该点距 EF< ， 则F为顶点。再找FH区段内距FH最远点位G，G距FH> ，贝U G为锚点；

⑧ 以FG、GH代替FH，先找GH区段内距FG最远点位，发现该点距 FH< ， 则G为顶点，再找GH区段内距GH最远点位，发现该点距 GH< ，则H为最 后一个顶点。

可见，算法每次迭代所傲的工作都是相同的： 找出拟合矢量区段内的最远点 位，求出该点到拟合矢量的距离，将这个距离与 比较，然后分两种情况，这个距

离若小于，则拟合矢量的末点为顶点，若大于 ，则为锚点，因此易于编程。 在寻找拟合矢量区段内的最远点位时，仅对该区间内的所有点位与拟合矢量之 间进行水平（或者垂直）方向的距离比较，不必求出每点到拟合矢量的真正距离。 4.5.2拟合函数

在得到特征点的基础上，我们建立了基于最小二乘法、Bezier曲线和B样条 的曲线

拟合模型。由于最小二乘法的原理较为通俗易懂，本文就贝奇尔曲线和B 样条的拟合方法作重点介绍。 •Bezier曲线拟合法

一条n次Bezier曲线［9］由n+1个控制点决定，依次连接这n+1个控制点即 得到特征多边形，如下图所示。

由Bezier曲线定义可得

v n

r n Jn,i u Vi

i 0

v

式中，n为Bezier曲线次数；

i为特征多边形顶点的序号，0 i n ；

u为参数，0 u 1;

Vi为特征多边形顶点的位置矢量； Jnj为伯恩斯坦基函数。

n i

Jn,i u end 1 u

其中，cn为组合数，

n! i! n i !

为表述简便，可将以上公式表示成矩阵的形式，即

V

n t,t,...,t,1 Mn C，…,Vn

nn 1

其中，Mn是n阶方阵，元素bij为

bj

1CnCni i

else

Bezier曲线的前提，由于已知 N+1个型值点

nijnij

已知控制点的位置矢量是使用

分别为Pj j 0,1,2,..., N，可用逆算法求解出控制点，即

n

J

V uu

n,i j i

u

V

p

i 0

写成矩阵的形式即为

V n 丄 n 1

r u t ,t ,..,t,1

1 C：Cn n 1 n 1 0

Cn Cn 1 10C°C°

n

1n1CnnCn

1 C n 1C1 1

Cn Cn 1

n 2

..... 0 ...

cc;

0 0

V V。 V V1 V Vn

该方程有n+1个未知数，N+1个方程，若n=N，贝昉程有唯一解；若n>N， 则方程有无穷多解，即控制点可随意选取；若 n根据上述分析可知，特征点的个数一般小于 得方程有唯一解。

边界条件：设与V u相邻的两条曲线如下，且u 0,1

4个，而且求得的Bezier曲线

要与两侧的曲线光滑拼接，故令n=N+2，将多余的两个方程用于边界条件的求解， 从而使

V V s u

m

J

m,i

u p

V

i 0

1

V w u

V

u Q

J|,i

i 0

k

①位置连续条件:

V V r 0 R V V

r 1 Q° ②斜率连续条件

V s v r , 1 A 1 V c r 0 v小 w 0

1时，即切矢共线，但模长不同，可满足数学

1时，拼接时切矢量相等，即曲率相等，此时拟合效果

其中，，为参数。当 上可到的条件；当 最好。 •B样条拟合法

与Bezier曲线拟合法类似，B样条曲线方程为

v

r u Ni,k u V

i 0

n

v

式中，N,k为B样条的基函数;

i为基函数序号; k

为基函数次数'

1 t

Ni,0 t

t 1

i

t

t

i 1

Ni,k t

N

t

0 else

ti k 1 t i,k 1

t

t

i k

t

i

Ni

1,k 1

i k 1

t

i 1

规定0

0。其逆运算与Bezier曲线类似,

Ni,k

v

Pj

同理，也多选两个控制点来补充端点条件

五、模型求解及结果分析

5.1模型求解

•Step 1:将示例位图导入MATLAB，并读取位图的像素灰度值;

•Step 2：在得到灰度值的基础上，选用合适的阈值确定方法，进行图像二值化； •Step 3:将二值化后的图像进行细化，得到单像素宽的图像轮廓骨架； •Step 4：将细化像素提取特征点及分割点；

•Step 5:用两个分割点的特征值进行 Bezire曲线拟合，得到有两个分割点的函 数方程； •Step 6:将所得到的分割点函数方程进行拼接，完成矢量化。 5.2模型结果及分析

实例位图如下图所示。

图示例原始位图

本例在求解过程中，使用了最大类间方差法和最大熵法两种方法来确定最佳 分割阈值，从而进行图像的二值化处理，两种方法得到的二值化图如下。

图 最大类间方差法二值化图 图 最大熵法二值化图

两种方法细化后的图像如下图。

图 最大类间方差法细化图

图 最大熵法细化图

为标记提取的间断点和分割点，我们将细化后的图进行颜色的黑白反转，使 背景为

白色，轮廓为黑色，以便更加明显地标记拟合需要的各类点。 颜色反转后 的图像如下。

由上图对比可知，两种方法得到的二值化图和细化图与原始图像的差距都不 大，且都较好的反映了原始图像的轮廓， 但是，最大类间方差法的二值化图像和 细化图像中的噪点比最大熵法的二值化图像和细化图像中的噪点多， 因此，为更 准确的进行特征点的提取和曲线拟合，我们采用最大熵法的图像进行后期处理。 提取间断点和分割点后的图像如下。

n Hi-4

图最大类间方差法细化颜色反转图 图最大熵法细化颜色反转图

图间断点分割点分布图

六、模型评价

6.1模型优点

1. 6.2模型缺点

1.

七、模型优化与推广 八、参考文献

[1] Pal S K ,Pal N R. Entropic thresholding [J].Signal Processing,19,16(2):97-108.

[2] 潘詰,吴一全.二维指数熵图像阈值选取方法及其快速算法[J].计算机应用, 2007,27⑷：982-985.

[3] Sahoo P K, Wilki ns C, Yeage J. Threshold selectio n using Re nyi's entropy[J].Pattern Recognition,1997, 30 ( 1) : 71-84.

[4] Abutaleb A S. Automatic Thresholdi ng of Gray-Level Picture Using Two-Dime nsional En tropies [J].Computer Visio n ‘Graphics ,and Image Processi ng ,19,47(1):22-32.

[5] Rao M ,Che n Y ,Vemuri B C,etal. Cumulative residual en tropy:A new measure of information [J].IEEE Trans.Inform.Theory,2004,50(6):1220-1228.

⑹郭臻，陈远知.图像阈值分割算法研究[J].中国传媒大学学报：自然科学版, 2008, 15(2): 77-82.

[7]

计算机图象处理中点阵图矢量化的研究及其在图形识别中的应用 微型机与应用,1988, 5: 004.

宁飞.[J].

[8] 王晓丹，吴崇明.基于MATLAB的系统分析与设计：图像处理[M].西安电 子科技大学出版社，2000.

[9] 计算机辅助几何造型技术[M].科学出版社，2009.