如何描述物体边界：常用方法介绍

区域外部形状是指构成区域边界的像素集合。为了描述目标物体的二维形状，除了前节所述的区域描述方法之外，通常采用的另一种方法是利用目标物的边界来表示物体，即所谓的边界描述。本节介绍几种常用的边界描述方法，包括边界的曲率、挠率、链码、傅立叶形状描述子、基于内角的形状特征和骨架化等。

（一）边界的曲率和挠率

如图4-24所示，当一点沿曲线以单位速率行进时，切向量t转动的快慢反映了曲线的弯曲程度，是用它对弧长的一阶导数衡量的。若设Δθ为切向量t（s）和t（s+Δs）间的夹角，s为弧长参数，那么有

在微分几何中，用曲率κ来刻画曲线在一点的弯曲程度，它等于曲线的切矢量相对

于弧长的转动率

设曲线方程为r=r（t），那么曲率的计算公式为

将上式展开，用r′和r″的分量表示，则为

κ随点在曲线上的移动而变化，如果κ≡0，那么该曲线是一条直线。

曲线的曲率是由曲线本身的形状决定的，与它的参数表达式和它在空间中的位置是无关的。反过来，这个量也可以唯一地确定平面曲线的形状。同时某个边界点曲率的正负也描述了边界在该点的凹凸性。

由于从数字图像中提取出的边界一般是离散且粗糙不平的，若直接计算边界的曲率，则结果不可靠，一般是用线段逼近边界后计算线段交点处的曲率，或者对边界进行曲线拟合，再根据式（4-82）和式（4-83）进行计算。

图4-24　切矢量沿弧长的变化

图4-25　副法向量沿弧长的变化

挠率用来表示空间曲线的扭曲程度。曲线在一点的切向量t和主法向量n所张的平面是密切平面，该平面的法向量是曲线的副法向量b，如图4-25所示。反映了密切平面方向转动的快慢，因而它刻画了曲线偏离平面曲线的程度，即曲线扭曲的程度。曲线在一点的挠率τ等于副法向量b（或密切平面）对弧长的转动率：

式中，Δθ是副法向量b（s）和b（s+Δs）之间的夹角，s是弧长参数。对于一般参数曲线r=r（t），挠率的计算公式为

将上式右端展开，用r′、r〃和r‴的分量表示为

挠率的符号按照如下规定选取：当一点沿曲线正向（即参数s增大的方向）移动时，b与n反向，则τ取正号；反之，则取负号。对于平面曲线，密切平面与曲线所在平面一致，因而副法向量是固定不变的，即b′≡0，故挠率τ≡0。

（二）链码描述符

区域边界曲线可用一组被称为链码的代码来表示，即Freeman方向链码。这种链码组合既利于有关形状特征的计算，也利于节省存储空间。

该方法采用曲线起始点的坐标和斜率（方向）来表示曲线。对于离散的数字图像而言，区域的边界轮廓可看作由相邻边界像素之间的单元连线逐段相连而成。对于某像素，把该像素和其8邻域内各像素的连线方向按八链码进行编码，即用0、1、2、3、4、5、6、7表示8个方向，其中偶数码为水平或垂直方向的链码，码长为1；奇数码为对角线方向的链码，码长为，如图4-26a所示。边界链码具有行进的方向性，如图4-26b所示，若以s为起始点，按逆时针的方向编码，所构成的链码为556570700122333，图4-26c的链码为65643432176010，若按顺时针方向编码，则得到的链码与逆时针方向的编码不同。用边界链码存储一个物体的边界，只需要一个起始点的（x，y）坐标以及每个边界点的三比特信息（8方向）或二比特信息（4方向）。

图4-26　八链码示意图

a）按顺时针方向的8个方向码　b）八链码示例1　c）八链码示例2

使用链码时，起点的选择是很关键的。对同一个边界，如果用不同的边界点作为链码起点，得到的链码是不同的。解决此问题的方法是：将这些方向数依1个方向循环，以使它们所构成的自然数的值最小。将如此转换后所对应的链码起点作为这个边界的归一化链码的起点，如图4-27所示。

图4-27　边界链码的起点归一化

a）原链码10103322　b）起点归一化链码01033221

根据链码可以计算一系列的形状特征，例如：

1.区域边界的周长

假设区域的边界链码为a1，a2，…，an，每个码段ai所表示的线段长度为Δli，那么该区域边界的周长为

式中，ne为链码序列中偶数码个数；n为链码序列中码的总个数。

2.区域的面积

按顺时针方向编码，边界链码对x轴的积分就是边界曲线所包围区域的面积S；

式中，yi=yi-1+ai2；y0是初始点的纵坐标；ai0和ai2分别是链码第i环的长度在k=0（水平）和k=2（垂直）方向的分量。对于封闭链码（即初始点与终点坐标相同），y0能任意选择。

3.对x轴的一阶矩

4.对x轴的二阶矩

5.形心（xc，yc）

式中是链码关于y轴的一阶矩。它的计算过程为：先将链码的每个方向码旋转90°，得到

然后利用式（4-89）进行计算。

6.两点之间的距离

如果链码中任意两个离散点之间的码为a1，a2，…，am，那么这两点之间的距离是

7.微分链码（差分码）

用链码表示目标边界时，如果目标平移则链码不会发生变化，而目标旋转则其链码会发生变化。为解决这个问题，可利用链码的一阶差分来重新构造一个表示原链码各段之间方向变化后的新序列，即差分链码，相当于把原链码进行了旋转归一化操作，因此差分链码具有旋转不变特性。如图4-28所示，左图的目标在逆时针旋转90°后成为右图的形状，原链码发生了变化，但微分链码并没有变化。左图中目标的原链码是10103322，微分链码是33133030；右图中目标的原链码是21210033，微分链码是33133030。

图4-28　链码的旋转归一化（利用一阶差分）

可以用两个相邻像素的码元方向数相减（后一个码元方向减去前一个码元方向），并对结果做模8运算得到微分链码。如图4-29所示，微分链码反映了边界的曲率，峰值处显示了凹凸性。图4-30为一个封闭曲线的链码和微分链码示例。

图4-29　链码和它的导数

a）目标区域及其边界　b）边界链码　c）微分链码

图4-30　封闭曲线的链码和微分链码

a）链码为5565707001223324，微分链码为017217101101072　b）链码为7707121223445546，微分链码为017217101101072

8.形状数（归一化差分码）

形状数是基于链码的一种边界形状描述符。根据链码的起点位置不同，一个用链码表达的边界可以有多个一阶差分。一个边界的形状数是这些差分中值最小的一个序列。例如图4-31中归一化前图形的基于4-方向的链码为10103322，微分码为33133030，形状数为03033133。

（三）傅立叶描述子

傅立叶描述子是区域边界变换的一种经典方法，在二维和三维形状分析中起着重要作用，其基本思想是用物体边界的傅立叶变换描述其形状，用较少的参数描述很复杂的边界。下面介绍两种计算物体边界的傅立叶描述子的方法。

1.方法1

可以用简单曲线来表示区域边界，设封闭曲线r在直角坐标系中表示为y=f（x），若对y=f（x）直接进行傅立叶变换，则变换的结果依赖于坐标x和y的值，不能满足平移和旋转不变性的要求。为了解决上述问题，引入以封闭曲线弧长为自变量的参数表示形式：z（l）=（x（l），y（l）），若封闭曲线的全长为L，则L≥l≥0。如图4-31所示，若曲线的起始点L=0，θ（l）表示曲线上某点的切线方向。设φ（l）为从起始点到弧长为l的点的切线的旋转角度，φ（l）随弧长l的变化而变化，显然它是平移和旋转不变的。则φ（l）=θ（l）-θ（0）。把φ（l）化为[0，2π]区间上的周期函数，用傅立叶级数展开，那么变换后的系数可用来描述区域边界的形状特征。因此，φ（l）的变化规律可以用来描述封闭曲线r的形状。

曲线r可看作多边形折线的逼近，假设折线有m个顶点ν0，ν1，ν2，…，vm-1，且该多边形边长νi-1νi的长度为Δli（i=1，2，…，m），则它的周长为

引入新的变量t；令，则γ∈[0，L]，t∈[0，2π]，定义

那么φ∗（t）为[0，2π]上的周期函数，且φ∗（0）=φ∗（2π）=0。φ∗（t）在封闭曲线r平移和旋转时均不变，并且φ∗（t）与r是一一对应的关系。由于φ∗（t）为周期函数，可用傅立叶系数对它进行描述，在[0，2π]上展开成傅立叶级数为

各谐波的系数分别为

图4-31　傅立叶描述子图解

2.方法2

把一个闭合边界曲线放到复平面上去，形成一个复数序列，即横坐标为实轴，纵坐标为虚轴。对该复数序列进行离散傅立叶变换，就得到轮廓的傅立叶描述。可以采用快速傅立叶变换来提高算法效率。

假设一个二维物体的轮廓是由一系列坐标为（xs，ys）的像素组成，其中0≤s≤N-1，N是轮廓上像素的总数。从这些边界点的坐标中可以推导出三种形状表达：

复坐标函数是用复数表示的边界像素坐标：

式中，（xc，yc）是物体的质心。对复坐标函数的傅立叶变换会产生一系列复数系数。这些系数在频率上表示了物体的形状，其中低频分量表示形状的宏观属性，高频分量表达了形状的细节特征。可以从这些变换参数中得出形状描述符。为了保持旋转无关性，仅保留参数的大小信息，而略去相位信息。通过将参数的大小除以直流分量的幅值，可保证缩放的无关性。(www.xing528.com)

轮廓线的曲率函数κ（s）可以表示为

式中，θ（s）是轮廓线的切向角度：

质心（重心）距离定义为从物体边界点到物体质心（xc，yc）的距离：

对于曲率函数和质心距离函数，由于其傅立叶变换是对称的，因而只考虑正频率轴。基于曲率函数的形状描述符表示为

式中，Fi表示傅立叶变换参数的第i个分量。类似的，由质心距离所导出的形状描述符为

对于复坐标函数，同时采用正频率分量和负频率分量。由于直流参数与形状所处的位置有关而不被采用。因此用第一个非零的频率分量来对其他变换参数进行归一化。复坐标函数所导出的形状描述符为

（四）基于内角的形状特征

首先，将物体的边界近似的表达成多边形的形式。多边形的内角

对形状的表达和识别非常重要。显然，基于内角的形状描述与形状所在位置、旋转和大小无关。以图4-32中内角θ=∠abc的计算为例，设a、b、c三点的中心为p，则有

式中，o为坐标原点。如p在多边形外部，则θ＞180°（图4-32a），此时

如p在多边形内部，则0°≤θ≤180°（图4-32b），此时

图4-32　物体边界的内角

a）凹形边界　b）凸形边界

以下是从内角导出的一系列形状特征：

1.顶点数

多边形的顶点数目越多，形状就越复杂。在识别图像中的目标物时，把具有不同顶点数目的两个形状当作不相似的两个形状是有一定合理性的。

2.内角平均值

多边形所有内角的平均值从一定程度上反映了多边形的形状属性。例如三角形的内角平均值为60°，与矩形的内角平均值90°之间有较大差别。

3.内角标准方差

多边形内角的标准方差为

式中，是内角的平均值。标准方差δ是多边形形状的总体描述：多边形越规则，δ值越小；反之，则δ值越大。因此，可以用值来分辨正多边形和不规则多边形。

4.内角直方图

首先将0°～360°的角度范围等分成k个区间，作为直方图的k个bin，然后统计每个角度区间中的内角数目，得到内角直方图，反映了内角的总体分布。

（五）骨架化

骨架是描述图像几何及拓扑性质的重要特征之一，骨架化是一种将区域结构形状简化为图形的重要方法。一个好的骨架应满足以下要求：

1）拓扑等价：骨架应和原始图像拓扑等价，即具有相同数量的前景目标、背景目标和孔。

2）细：骨架为单像素宽。

3）居中：骨架应位于目标区域的中心。

1.距离变换

自Rosenfeld和Pfaltz于1966年首次提出距离变换的概念以来，它已被广泛应用于图像分析、计算机视觉和模式识别领域中，一般用来实现目标细化、骨架抽取、形状的插值和匹配、粘连物体的分离等。

距离变换是针对二值图像的一种变换。二值图像可以认为仅包含目标和背景两类像素，目标像素的灰度值为1（即255），背景像素的灰度值为0。如图4-33所示，设P为灰度值为1的像素区域，Q为灰度值为0的像素区域，求从P中任意像素到Q的最小距离的处理叫作二值图像的距离变换。距离变换即求二值图像中各1像素到0像素的最短距离的处理，是把任意图形转换成线划图的最有效方法之一，其结果不是另一幅二值图像，而是一幅灰度级图像，即距离图像，图像中每个像素的灰度值为该像素与距其最近的背景像素之间的距离。

图4-33　二值图像中1像素和0像素之间的距离

现有的距离变换算法主要采用两类距离测度来度量任意两个像素之间的距离：非欧氏距离和欧氏距离。前者常用的有街区距离、棋盘距离、倒角距离等，算法采用串行扫描实现距离变换，在扫描过程中传递最短距离信息。这些算法简单快速、易于实现，但得到的仅是欧氏距离变换（Euclid Distance Transform，EDT）的一种近似值，在很多应用中不能满足精度要求，必须使用欧氏距离变换。

二维图像平面上的两点（x1，y1）和（x2，y2）之间的欧氏距离表示为

在二值图像中，1代表目标点，0代表背景；在灰度图像中，像素的灰度值表示该像素到最近目标点的距离值。这样一幅大小为M×N像素的图像可以表示为一个二维数组A[M][N]，其中A[i][j]=1对应的像素表示目标点，A[i][j]=0对应的像素表示背景点。设B={（x，y）|A[i][j]=1}为目标点集合，则欧氏距离变换就是对A中所有的像素点求最短欧氏距离：

式中，Distance[（i，j），（x，y）]为点（i，j）和（x，y）之间的欧氏距离，从而得到二值图像A的欧氏距离变换图。

以4邻接方式为例，对于二值图像f（i，j），将其距离变换k次的图像为gk（i，j），

对全部i和j取gk+1（i，j）=gk（i，j）时，gk便是所求的距离变换图像。当f（i，j）=1时，g0（i，j）是一个非常大的值；f（i，j）=0时，g0（i，j）=0。

在经过距离变换得到的图像中，最大值点的集合就形成骨架，即位于图像中心部分的像素集合，也可以看作是图形各内接圆心的集合，它反映了原图形的形状，如图4-34所示。给定距离和骨架就能恢复该图形，但恢复的图形不能保证原始图形的连接性。该方法常用于图形压缩、提取图形幅宽和形状特征等。

图4-34　对一幅二值图像进行距离变换

a）原始二值图像示意图　b）距离变换结果

2.中轴变换

中轴变换（Medial Axis Transform，MAT）是一种用来确定物体骨架的细化技术，在不影响原图拓扑性的基础上，通过抽取表达原图形状的最关键点，使得原图中宽度大于1个像素的线条变成单像素，每个骨架点保持与边界点距离最小。物体内部的一点位于其中轴上的充要条件是，它是与物体边界相切于两个相邻点的圆的圆心，如图4-35所示。

中轴变换的算法和第3章中介绍的细化算法有些相似，也是采取逐次去除边界的方法来进行的，同时中轴变换也不能破坏图像的连通性。具有边界B的区域R的MAT是按如下方法确定的：如图4-36所示，对于R中的点p，在B中搜寻与它最近的点，即找p在B上“最近”的邻居。如果对p能找到多于一个这样的点（即有两个或两个以上的B中的点与p同时最近），就可认为p属于R的中轴线或骨架，或者说p是一个骨架点。

图4-35　物体的中轴示意图

图4-36　中轴变换示意图

目前，中轴变换的方法很多，通常利用二值形态学操作判断像素的8邻域情况。例如，可通过对二值图像进行两次删除，每次删除满足4个条件，完成中轴变换。即设置一个5×5的模板S，S中各元素的取值取决于模板所对应图像中的不同像素，如果S某一个位置所对应的像素值为255，将模板上该位置赋为0，否则赋为1。设N（S[2][2]）表示以S[2][2]为中心的3×3邻域内目标像素（即黑点）的个数，以S[2][2]为中心取3×3邻域，则T（S[2][2]）表示序列{S[1][2]，S[1][1]，S[2][1]，S[3][1]，S[3][2]，S[3][3]，S[2][3]，S[1][3]，S[1][2]}中由0变1的次数。以S[1][2]为中心取3×3邻域，则T（S[1][2]）表示序列{S[0][2]，S[0][1]，S[1][1]，S[2][1]，S[2][2]，S[2][3]，S[1][3]，S[0][3]，S[0][2]}中由0变1的次数。以S[2][1]为中心取3×3邻域，则T（S[2][1]）表示序列{S[1][1]，S[1][0]，S[2][0]，S[3][0]，S[3][1]，S[3][2]，S[2][2]，S[1][2]，S[1][1]}中由0变1的次数。第一次删除的4个条件是：条件1：2≤N（S[2][2]）≤6；条件2：T（S[2][2]）=1；条件3：S[1][2]×5[2][1]×S[3][2]=0；条件4：S[2][1]×S[2][3]×S[3][2]=0。第二次删除的4个条件中前三个条件与第一次删除相同，将条件4改为：S[1][2]×S[2][3]×S[3][2]=0。先判断第一组的4个条件，遍历整幅图像，如果有像素满足全部条件则删除该点，然后再用处理后的图像进行第二次判断，遍历整幅图像，如果有像素满足第二组的4个条件，则删除该点。对处理后的图像重复进行上面的两步操作，直至没有点可以删除，这时剩下的点就是图像的中轴。

算法的具体步骤如下：

1）获取原图像的首地址以及图像的高度和宽度。

2）开辟一块内存缓冲区，并初始化为255。

3）对图像进行遍历，如果当前像素的灰度值为255（即背景），则跳过该像素。

4）如果当前像素的灰度值为0（即物体），则定义一个5×5的结构元素，计算结构元素中各个位置上的值，为防越界，不处理外围的2行、2列像素，从第3行第3列开始判断，将结构元素中心与当前像素重合，如果结构元素所覆盖的位置下，像素的灰度值为255，则在结构元素上同样的位置处置0，否则是目标，置1。

5）依次判断第一组的4个条件，若4个条件全部满足则删除该点，否则判断下一个像素，直至全部像素处理过一遍。

6）对处理过的图像依次判断第二组的4个条件，若4个条件全部满足则删除该点，否则判断下一个像素，直至全部像素处理过一遍。

7）循环执行步骤5）和6）直至没有像素可以删除。

8）将结果保存到内存缓冲区。

9）将结果由内存缓冲区复制到原图的数据区。

在进行图像识别时，首先对被处理的图像进行中轴变换有助于突出形状特征和减少信息量的冗余。中轴变换利于找出细长而弯曲物体的中心轴线，例如图4-37是用欧氏距离计算出的一些区域的骨架，图4-38是对血管骨架的提取。其他的形状描述子，如物体具有的分支数和物体的总长，可以从中轴变换图计算出来。对二值图像来说，中轴变换能够保持物体的原本形状，因此该变换是可逆的，物体可以由它的中轴变换而得到重建。

图4-37　用欧氏距离计算出的一些区域的骨架

图4-38　X射线冠状动脉造影图像中主要血管分支骨架的提取示例