基于像素-亚像素级形态分析的路面三维图像裂缝自动识别算法

彭 博, 黄大荣,郭 黎, 蔡晓禹, 李少博

(1. 重庆交通大学 交通运输学院, 重庆 400074; 2. 湖北民族学院 信息工程学院, 湖北 恩施 445000)

0 引 言

裂缝是路面的主要破损形式之一，它严重影响路面质量、服务性能与使用寿命，需要对其进行及时有效的检测分析，以评估路面性能、辅助养护管理决策、改善路面结构和材料设计。自1988年中国建成第一条高速公路以来，我国高等级公路经历了近30年的快速建设期，至今相关规划建设工作依旧如火如荼地展开。在新建、改扩建公路特别是高等级公路的同时，路面养护管理的紧迫性和重要性也日渐凸显出来。因此，路面裂缝自动检测与分析已然成为路面工程研究领域和道路养护管理部门的关注焦点。

为此，国内外学者对路面破损检测技术进行了长期研究，如效率低、劳动强度高的人眼观测法、录像测度法，能检测路面厚度、裂缝、空洞的探地雷达法[1]，能全面地采集路面裂缝、松散、坑槽等信息且具有直观观察优势的摄像测量法[2]。其中，摄像测量法得到了长足的发展，自20世纪70年代以来，出现许多基于摄像测量技术的路面破损检测系统，如日本的Komatsu系统[3]、瑞典的PAVUE系统[4]、美国的PCES系统[4]和DHDV系统[5]、中国的CiCs系统等。从采集的路面原始数据而言，大多检测系统可获取路面高清灰度图像(如PCES系统、CiCs系统)，随着激光摄像技术的发展，近年来少部分系统可获得路面高精三维深度图像(如第三代DHDV系统可获取1 mm精度的路面三维图像)。

对于路面裂缝识别，在路面高精三维图像出现以前，研究人员基于路面灰度图像提出了许多识别算法，如阈值分割算法[6-7]、边缘识别算法[8-9]、种子识别算法[10-11]、频域变换算法[7,12]、纹理分析方法[13]、机器学习算法[14]等。由于灰度图像容易受到车道标线、阴影及油污等影响，会降低路面图像质量、增加裂缝检测难度，而路面三维深度图像真实体现了路面深度信息与形态轮廓，可有效规避灰度图像的缺陷。于是，随着路面高精三维图像的出现，部分人员开始研究基于三维图像的裂缝自动识别算法，提出了种子生长法[15-16]、各向异性测度算法[17]、匹配滤波法[18]、光影模型法[19]等，一些既有的灰度图像处理理论与方法也被融合到了三维图像裂缝识别过程中。

然而，尽管算法众多，由于路面状况的复杂多样性，各个算法通常具有较苛刻的适用条件，无法普遍适用于各类复杂路面的裂缝识别，准确快速识别完整的路面裂缝图像仍然是一项挑战。

为了进一步推动路面裂缝的高精快速自动识别，区别于既有算法利用像素级或子块级(q像素×q像素区域视为一个子块)图像信息提取裂缝，笔者基于像素-亚像素级图像信息融合和图像骨架-区域形态分析，提出了路面高精三维图像的裂缝自动识别算法。然后，编程实现该算法。最后，基于同一路面图像集，对本文算法和另外5种既有算法进行了性能测试分析，探讨了算法适用性。

1 算法简介

本文算法以融合边缘检测和区域生长图像为基础，基于候选目标图像的骨架形态特征和连通域形态特征进行裂缝自动识别。数据来源于DHDV系统和国内某公司采集的1 mm/像素路面三维图像(992像素×992像素，像素值代表深度)，如图1。算法开发和测试环境为Microsoft Visual Studio 2010+OpenCV2.4.9+Halcon11.0。

图1 1 mm/像素路面三维图像Fig. 1 3D pavement image of 1mm/pixel

算法基本流程如图2。

图2 方法总流程Fig. 2 Overall procedures of the proposed method

算法主要步骤如下：

1)图像分割。分别应用边缘检测和区域生长分割方法提取候选目标图像，并进行融合处理，获得分割图像。

2)图像骨架分析与处理。提取分割图像的亚像素级骨架图像和像素级骨架图像，然后分别展开图像骨架重构分析。该步骤生成重构的骨架图像，可增强候选目标图像的连续性和完整性。

3)基于形态分析的裂缝自动识别。利用形态学操作和连通域轮廓长度、面积、圆度、外接圆半径、等效椭圆轴长等形态特征，融合骨架图像、识别路面裂缝。

2 图像分割

2.1 Canny边缘检测

对于图像目标识别而言，边缘检测可大幅度地减少数据量，剔除大量不相关或弱相关信息，因此，许多研究均以边缘特征进行目标识别。

1986年，J. CANNY[8]提出了一种至今仍广泛应用的经典边缘检测算法，该算法主要思想和步骤为：①应用高斯滤波消除图像噪声；②通过Sobel算子卷积获得图像梯度强度与方向；③应用非极大抑制技术(Non-maximum suppression)保留每个像素点上梯度强度的极大值，而删掉其他的值；④通过双阈值技术界定弱边界和强边界，并进行边缘连接处理、形成边缘图像。

对图1进行Canny边缘检测，获得边缘强度图像；然后利用简单阈值分割法或OTSU分割法[6]，得到二值边缘图像E1，如图3。

图3 边缘检测图E1Fig. 3 Edge detection image, E1

2.2 区域生长分割

图像分割是数字图像处理中的关键技术之一，它是根据图像的灰度、空间纹理、几何形状等特征把图像划分成若干个不互相交的区域，使得这些特征在同一区域内，表现出一致性或相似性，而在不同区域之间表现出明显的不同。

区域生长[20]是一种常用的图像分割方法，其基本思想是将具有相似性质(如灰度值)的像素集合起来构成区域，区域生长的首要任务是选择有意义的特征和确定相似性准则。一般可以像素灰度值为特征，以灰度值之差衡量相似性。

以任意两个相邻n像素×n像素矩形为例，令两个矩形的特征值(可以是中心像素灰度值，亦或n×n个像素灰度均值)分别为g1和g2，若满足以下条件，则两个矩形合并为一个区域：

|g1-g2|

(1)

式中:T为阈值参数。

对图1进行区域生长分割，得到区域生长图E2，见图4。

图4 区域生长图E2Fig. 4 Region growing image, E2

2.3 图像融合

Canny边缘检测与区域生长基于不同的思想获得了候选目标图像，为了充分利用它们提取的目标信息、增强目标图像的完整性，将E1、E2进行融合处理(即取并集)，得到融合分割图像G，如图5。

(2)

图5 融合分割图GFig. 5 Fused segmentation image, G

3 图像骨架分析与处理

图像骨架指图像目标的中轴，可以有效反映出原物体形状的连通性和拓扑结构。本文算法基于图像分割结果，提取和重构图像骨架(包括像素级和亚像素级骨架图像)，为裂缝识别奠定基础。

3.1 像素级图像骨架分析

基于融合分割图像G，利用图像细化算法[21]获取像素级骨架图像G1。然后分两步进行膨胀处理，其目的是增强图像的连续性和完整性：①对G1进行膨胀处理，并再次提取像素级骨架图像G2；②对G2中所有长度大于L1(笔者取L1=40像素)的骨架曲线进行膨胀处理，得到图像G3，如图6。

图6 像素级图像骨架处理结果G3Fig. 6 Processing result of pixel level image skeleton, G3

3.2 亚像素级图像骨架分析

3.2.1 亚像素级图像骨架生成

亚像素(Subpixel)[22]的概念如图7，每个物理像素可以分为m×m个亚像素，用数学插值算法可计算每个亚像素值，通过亚像素图像处理可更大限度利用图像信息。

图7 亚像素Fig. 7 Subpixel

基于宽度为1像素的像素级骨架图像G1，每个像素的中心亚像素相连即可构成亚像素骨架G4。以G1的局部图像为例，亚像素骨架生成如图8。

图8 亚像素骨架生成Fig. 8 Generation of subpixel skeleton

3.2.2 亚像素级图像骨架重构

对于亚像素级骨架图像G4的骨架曲线，分两步进行连接处理。①对于任意相邻骨架曲线，若其端点距离小于阈值d1(笔者取d1=2像素)，则连接这两个端点；②对于所有长度大于L2(笔者取L2=20像素)的骨架曲线，若其端点距离小于阈值d2(笔者取d2=10像素)，则连接这两个端点。

以图8亚像素骨架图像为例，连接结果见图9。由图9可以看出，骨架连接处理能有效增强目标图像结构的封闭性、连续性和完整性，有助于裂缝识别。由于骨架连接结果仍为亚像素级图像，为便于开展后续识别工作并降低运算成本，还需将其还原为像素级图像，如图9。

图9 骨架重构过程Fig. 9 Skeleton reconstruction process

通过步骤①、②处理后进一步生成像素级骨架图像G5，如图10。

图10 骨架重构图像G5Fig. 10 Reconstructed skeleton image, G5

4 基于形态分析的裂缝自动识别

主要利用形态学运算(膨胀和闭运算)，同时结合连通域的轮廓长度、面积、圆度、最小外接圆半径、等效椭圆长短轴、扁平率等多种形态特征，对骨架图像进行融合处理，提取裂缝目标。

4.1 形态学算子简介——膨胀和闭运算

膨胀是利用t×t的结构元素，将与前景物体接触的所有背景点合并到该物体中，使边界向外部扩张，可以用来填补物体中的部分空隙。腐蚀则与膨胀相反，是利用t×t的结构元素使边界向内部收缩，可以用来消除小且无意义的物体。先膨胀后腐蚀的过程称为闭运算，可用来填补物体内细小空洞、连接邻近物体、平滑其边界。如图11，左图为右图的腐蚀效果，反之，右图为左图的膨胀效果。闭运算效果如图12。

图11 腐蚀和膨胀Fig. 11 Erosion and dilation

图12 闭运算Fig. 12 Closed operation

4.2 连通域形态特征

在图像中，最小的单位是像素，每个像素周围有8个邻接像素，常见的邻接关系有2种：四邻接与八邻接。四邻接一共4个点，即上下左右；八邻接的点一共有8个，包括了对角线位置的点，如图13。

图13 四邻接与八邻接Fig. 13 4-neighbourhood and 8-neighbourhood

根据四邻接或八邻接关系，彼此连通的像素点形成的区域称为连通域。在图像处理中，常常对连通域展开计算分析，以达到图像分割、噪声去除和目标识别等目的。笔者基于八邻接关系建立连通域，并展开长度、面积、圆度和等效椭圆长短轴等形态特征计算分析。

1)轮廓长度L：连通域轮廓边缘所含像素总数。

2)面积S：连通域面域所含像素总数。

3)圆度R[23]：接近理论圆的程度，连通域棱角越多越尖锐则圆度越差，反之棱角圆滑，圆度就高。

4)最小外接圆半径r。

5)等效椭圆[24]长短轴：每个连通域可求得一个对应的椭圆，其方程与连通域的转动惯量形式一致。假定这个椭圆区域是均质的且椭圆与连通域的转动惯量相同，那么椭圆的参数就能反映区域的特性，这个椭圆称为等效椭圆。据此可求得等效椭圆长轴ra和短轴rb。见图14。

图14 连通域等效椭圆Fig. 14 Equivalent ellipse of connected component

6)扁平率w：裂缝一般呈现为狭长的带状，定义扁平率计算公式如式(3)。可见，w越大，连通域越狭长。

(3)

4.3 图像骨架融合

对像素级图像骨架分析结果G3和亚像素级图像骨架重构图像G5进行融合，其主要思想是通过设置轮廓长度、面积、最小外接圆半径和圆度阈值，保留轮廓长度、面积和最小外接圆半径较大且圆度较小的曲线对象。通过合理设置阈值，可保留大型目标、剔除尺寸明显较小的曲线对象，避免误删较小的路面裂缝。

融合结果为G6，如图15。

图15 骨架融合图像G6Fig. 15 Fused skeleton image, G6

融合处理算法伪代码如下。

算法1 图像骨架融合处理算法(“←”表示赋值)

1：intN1,N2←骨架图像G3和G5的曲线数量

2：FOR inti=1 toN1

3：G3[i] ←G3第i个骨架曲线

4：L3[i] ←G3[i]的轮廓长度

5：IFL3[i]

6：DELETEL3[i] // 删除曲线L3[i]

7：END IF

8：END FOR // 得到骨架图像G3’

9：FOR inti=1 toN2

10：G5[i] ←G5第i个骨架曲线

11：L5[i],S5[i]←G5[i]的轮廓长度，面积

12：R5[i],r5[i]←G5[i]的圆度，最小外接圆半径

13：IFL5[i]R0 ORr5[i]

14：DELETEL5[i]// 删除曲线L5[i]

15：END IF

16：END FOR//得到骨架图像G5’

17：G6←∪{G3′,G5′}//取G3′和G5′并集，得到融合图像G6

4.4 裂缝识别

裂缝识别主要包括2个步骤：①对骨架融合图像G6进行闭运算和膨胀操作，获得连通域图像如图16。②应用轮廓长度、面积、圆度等形态学特征筛选满足要求的连通域，最终裂缝识别结果见图17。

图16 连通域图像Fig. 16 Connected component image

图17 裂缝二值图Fig. 17 Binary crack image

裂缝识别算法伪代码如下。

算法2 基于形态分析的裂缝识别算法(“←”表示赋值)

1：Image0←G6闭运算和膨胀操作结果，见图16

2：M←Image0的连通域个数

3：C[i]←Image0的第i个连通域

4：L[i],S[i] ←C[i]的轮廓长度，面积

5：R[i],r[i] ←C[i]的圆度，最小外接圆半径

6：ra[i],w[i] ←C[i]的等效椭圆的长轴长度，扁平率

7：FOR inti=1 toM

8：IFL[i]

9：DELETEC[i] //移除连通域C[i]

10：END IF

11： IFR[i]>R1 ORr[i]

12：DELETEC[i] // 移除连通域C[i]

13：END IF

14： IFra[i]

15：DELETEC[i] // 移除连通域C[i]

16：END IF

17：END FOR // 获得裂缝识别结果二值图，如图17

5 算法性能测试

5.1 测试方法

为了客观定量地评价算法性能，首先，基于150张1 mm/像素路面三维图像(992像素×992像素)手动标注出真实裂缝图像，以此作为测试基准图像库。然后，分别应用各向异性测度算法[17]、Canny边缘检测算法[8]、种子并行生长算法[16]、GAVILN算法[11]、OTSU分割方法[6]和本文算法对上述150张图像进行裂缝检测。最后，对比检测结果和人工绘制的真实裂缝图像，计算裂缝识别的准确率(Pr，Precision，反映检测结果的准确程度)、召回率(Re，Recall，反映检测结果的全面性)和F值，据此评价算法性能优劣[15-16]。

(4)

5.2 测试结果分析

图18展示了上述6个算法准确率、召回率和F值的累计频率分布曲线(150张图像的测试结果)。

图18 算法性能对比Fig. 18 Comparison of algorithm performance

从图18(a)可以看出，本文算法准确率较高，均值达90.45%，其余算法准确率由高到低依次为种子并行生长算法(均值44.34%)、各向异性测度算法(均值37.62%)、GAVILN算法(均值37.42%)、Canny检测(均值21.37%)和OTSU(均值3.95%)。

从图18(b)可以看出，本文算法召回率较高，均值为96.49%，其余算法召回率由高到低依次为Canny检测(均值56.41%)、种子并行生长算法(均值41.92%)、OTSU(均值40.30%)、GAVILN算法(均值38.75%)和各向异性测度算法(均值30.74%)。

F值综合考虑了准确率和召回率，其分布曲线如图18(c)，可以看出，F均值排序为：本文算法(均值90.72%)>种子并行生长算法(均值39.65%)>GAVILN算法(均值33.46%)>各向异性测度算法(均值30.32%)>Canny检测(均值25.85%)>OTSU(均值5.85%)。

对于本文算法的准确率和召回率，分布范围均为0～100%，其中，约有10%的样本准确率低于0.8、召回率低于0.9。换言之，绝大部分图像裂缝识别的准确率和召回率都能获得较高的水平。基于本文基准图像库对6种算法的测试结果显示，本文算法F值最高，表明它能准确而全面地识别路面裂缝。

图19为裂缝检测实例(宽992像素×高992像素的路面图像，实际覆盖范围约为1平方米)。OTSU分割方法将深度值较低的线性特征和区域均判为裂缝，因此召回率较高而准确率较低。Canny能检测深度阶跃变化的边缘特征并将其视为裂缝，召回率较高，但噪声干扰较多，准确率下降。通常而言，OTSU和Canny可作为裂缝识别算法中的一个环节。种子并行生长算法和GAVILN算法都是基于种子生长的裂缝识别方法，将原始路面图像划分为若干个y像素×y像素的正方形子块(y一般取4或8，每个小方块将被识别为裂缝种子或非裂缝子块)，然后进行裂缝种子识别与连接，因此，这两种算法识别的裂缝最小单位为子块级，容易漏检细小的裂缝，如图18(d)、图18(e)第3～6张图；此外，由于种子并行算法具有并行融合机制，检测结果总体上优于GAVILN算法。各向异性测度算法融入了像素级与子块级裂缝检测信息，不过基于图像子块进行了去噪处理，导致许多细小裂缝和裂缝片段被滤除，因此漏检情况较严重。

图19 裂缝识别结果对比Fig. 19 Comparison of crack detection results

综合来看，本文算法之所以检测效果优于种子并行生长算法、GAVILN算法和各向异性测度算法，主要有两方面原因：①本文算法融合了像素级和亚像素级图像信息，适用于检测较细小的裂缝目标；而后三种算法融合了像素级与子块级图像信息，更适合于检测较粗较明显的裂缝目标；②本文测试图像含有较多细小裂缝。因此，算法评测的优劣较大程度上受到测试图像的影响，在开发裂缝自动识别算法的同时还需要分析其适用条件，明确适用范围。

6 结 语

笔者针对1 mm/像素的路面深度图像提出了基于像素-亚像素级形态分析的裂缝自动识别算法，主要包括了图像分割、像素级和亚像素级图像骨架分析、图像骨架融合与形态学分析处理等步骤。测试分析表明：本文算法融合了像素级和亚像素级图像信息，对于细小裂缝较多的路面图像，识别效果良好，优于经典的OTSU分割方法和Canny边缘检测，优于种子并行生长算法、GAVILN算法和各向异性测度算法。

裂缝自动识别是道路基础设施智慧养护与管理的重要基础，由于路面状况的复杂多样性，各种算法检测结果仍有准确率和召回率很低的情况，准确稳定地识别裂缝破损依然是一项挑战。由于裂缝识别算法均有特定的适用条件，因此，有必要对不同特性的路面图像开发针对性的算法，进而形成路面裂缝检测的算法库。同时，算法开发过程中，还需有效结合路面二维和三维特征，引入人工智能和深度学习理论方法，以提高裂缝识别的自动化与智慧化水平。