CN104809756B - 基于X-ray CT图像的沥青混合料空隙空间结构重构方法 - Google Patents

Abstract

本发明公开了一种基于X‑ray CT图像的沥青混合料空隙空间结构重构方法，为应用于虚拟数值分析的沥青混合料中空隙的细观结构提供一种基于X‑ray CT图像的沥青混合料空隙空间结构重构方法，实现混合料中空隙三维建模的自动化。本发明基于沥青混合料细观结构中的空隙分布统计，完成从沥青混合料试件X‑ray CT断面图像到空隙空间结构实体模型的自动转换，以辅助统计沥青混合料内部空隙的真实空间分布，并为基于真实细观结构的沥青混合料细观三维数值建模奠定基础。

Description

基于X-ray CT图像的沥青混合料空隙空间结构重构方法

技术领域

本发明涉及数值分析领域，具体是一种基于X-ray CT图像的沥青混合料空隙空间结构重构方法。

背景技术

沥青混合料由集料、沥青胶浆和空隙组成，其力学行为由三相组分材料的体积含量和空间分布共同决定。其中空隙对沥青路面的抗车辙能力、透水性、疲劳寿命、强度等各性能有着重要的影响。因此研究沥青混合料中空隙的含量及分布，对于准确预测道路的路用及安全性能十分关键。为了准确描述沥青混合料中空隙的形状属性和空间分布，以便在此基础上重构混合料细观数值模型并获得可信的分析结果，建立空隙空间结构的三维模型十分必要。在空隙结构三维建模中存在两个关键问题：一是选取何种参数作为空隙的建模参数；二是如何获取建模参数。

目前对沥青混合料中空隙的研究可以分为宏观和细观两个方面。宏观研究以空隙率为对象，分析实验成型条件对空隙率的影响规律以及空隙率对沥青路面路用性能的影响规律等。而单一空隙率不能反映空隙在混合料中的形状及分布属性，因此难以准确预测混合料的路用性能。细观研究通过图像处理技术分析混合料断面图像中的空隙分布，一方面运用统计方法获得空隙在混合料中的数量、分布及体积等三维参数，另一方面将空隙像素映射为有限元或离散元单元反映到混合料数值模型中模拟细观力学行为。其中，由二维信息统计获得的空隙三维参数由于并没有将空隙的空间分布结构作为参照，准确度有限；而将空隙像素映射为有限元或离散元单元而构建的混合料数学模型，也没有表现出空隙的真实形貌，难以准确模拟空隙对集料颗粒和沥青胶浆之间细观力学行为的影响。同时，其所运用的图像处理技术对空隙的识别率较低。

长安大学2008年国家发明专利《采用CT技术分析多孔沥青混合料空隙结构的方法》(专利号：ZL200810232287.5)利用X-ray CT设备获得混合料横向和纵向断面图像，通过图像处理技术分割出断面图像中的空隙像素后，统计空隙的数量、面积、周长、等效直径、空隙率等二维及三维参数。但是，所采用的基于 阈值分割的空隙像素提取方法识别精度较低，且依据二维分析结果而推导出的空隙三维参数可信度也不高。

华南理工大学2010届博士学位论文《基于X-ray CT和有限元方法的沥青混合料三维重构与数值实验研究》将沥青混合料X-ray CT扫描图像叠加获得混合料的三维可视化模型，设定灰度阈值区分三相材料组分像素后，将各像素映射为有限元单元，建立了混合料的有限元模型。但其构建的有限元模型中没有表现出空隙等各相组分的真实形貌，因此获得的细观力学仿真结果与试验数据偏差较大。

基于上述分析，目前对沥青混合料中空隙的研究方法，不能重现空隙真实的细观结构，而导致所建混合料数值模型的细观力学仿真结果与试验数据间的较大偏差。本发明基于图像处理与三维建模技术，完成对混合料中空隙空间结构的三维实体建模。

发明内容

本发明的目的是提供一种基于X-ray CT图像的沥青混合料空隙空间结构重构方法，应用于虚拟数值分析的沥青混合料中空隙的细观结构，实现混合料中空隙三维建模的自动化。

为了达到上述目的，本发明所采用的技术方案为：

基于X-ray CT图像的沥青混合料空隙空间结构重构方法，其特征在于：包括以下步骤：

(1)输入沥青混合料试件的X-ray CT断面扫描图像文件组graph[n]，用以重构空隙空间结构的实体模型；

所述沥青混合料试件是指以沥青材料和集料按适当比例配制，胶结成整体的圆柱形样本，记半径为radius，高度为height；

所述X-ray CT断面扫描图像文件组是指使用X-ray CT设备以一定纵向长度间隔扫描沥青混合料试件获取的图像数组，其中，图像数组中的元素均为BMP格式、分辨率为512×512像素的灰度图像文件；

所述灰度图像文件是指将彩色图像中像素的颜色分量值按照一定比例进行加权，并将加权值附于每个颜色分量而获得的图像文件，其中，灰度图像中所有像素的灰度均介于0到255之间；

所述空隙空间结构是指反映沥青混合料试件中所有空隙的几何形状及分布信息的空间结构；

所述实体模型是指以边界表示法表示的三维几何结构，记为模型M＝(FACE,EDGE,VERTEX,R)，其中FACE表示模型中面的集合，EDGE表示边的集合，VERTEX表示顶点的集合，R表示模型中各元素之间的邻接关系； 记录着面的类别：平面、样条面，以及面的曲面方程，其中face为FACE中的元素；

记录着边的曲线方程，其中edge为EDGE中的元素；

记录着顶点的坐标p(x,y,z)，其中vertex为VERTEX中的元素；

R＝{r₁,r₂}，其中：

r₁＝{(face₁,face₂,edge)|edge∈face₁∩face₂,edge∈EDGE,face_i∈FACE,i＝1,2}，r₁表示face₁和face₂相邻于边edge；

r₂＝{(edge₁,edge₂,vertex)|vertex∈edge₁∩edge₂,vertex∈VERTEX,edge_i∈EDGE,i＝1,2}，r₂表示edge₁和edge₂相邻于点vertex；

(2)使用空隙轮廓检测算法，识别出每幅X-ray CT断面扫描图像graph[i](i＝1,2,…,n)中各单个空隙的轮廓像素，获得相应的空隙轮廓面实体模型，并建立沥青混合料中各单个空隙与相应的轮廓面实体模型间的映射关系；

(3)输入各单个空隙所对应的全部轮廓面实体模型(k＝1,2,…,m)，调用ACIS7.0的底层建模功能Skinning，获得BODY类型的单个空隙实体模型M_k(k＝1,2,…,m)；

(4)输入各单个空隙实体模型{M_k|k＝1,2,…,m)}，调用ACIS7.0的底层建模功能Boolean_Unite，获得BODY类型的空隙空间结构实体模型M，该模型以边界表示法表示各单个空隙的三维几何结构和在沥青混合料试件中的分布规律，最终存储为SAT格式的图形文件，可应用于沥青混合料内部空隙的真实空间分布 统计及沥青混合料细观数值分析；

所述空隙空间结构实体模型M是指能够反映沥青混合料试件中一定尺寸之上的所有空隙的形状属性及在试件中分布规律的空隙空间结构三维模型；

所述沥青混合料细观数值分析是指基于沥青混合料细观结构、用于预测沥青混合料的力学和热传导等性能的数值仿真。

所述的基于X-ray CT图像的沥青混合料空隙空间结构重构方法，其特征在于：步骤(2)中，空隙轮廓检测算法描述如下：

1)对graph[i]进行纵横向间距均为3个像素的均匀方形网格划分，获得网格节点像素集合NS_i{node_p|p＝1,2,…,16384}，其中各网格节点像素按其在graph[i]中的位置自上而下、自左至右进行顺序编号，令p＝1、graph[i]中所有像素的已访问标志bV＝false；

2)使用空隙像素判别算法判断node_p是否为位于某一空隙轮廓内的空隙像素，若是，将node_p及其八邻域像素分别作为该空隙轮廓对应灰度差分八叉树 的根结点及根结点的子结点，置此9个像素的已访问标志bV＝true并置根结点的所有子结点的可扩展标志bE＝true后计算中每条连接线段的灰度差分值d_FC，令当前层号h＝2后执行步骤3)；若否，仅置node_p的已访问标志bV＝true，使p增1后执行步骤2)；

3)对第h层上的结点，依次判断每个结点与其父结点的连接线段对应的d_FC是否超过灰度差异阈值若是，则该结点对应的像素为空隙轮廓像素，将该像素存入空隙轮廓元中的轮廓像素链表并设置该结点的可扩展标志bE＝false；

4)对第h层上bE＝true的结点，将其对应像素的八邻域像素中bV＝false的像素作为该结点的子结点并置各子结点所对应像素的bV＝true，使h增1；

5)循环执行步骤3)～步骤4)，直至h等于预设最大层号H或第h层上的结点所对应像素均为空隙轮廓像素，完成构建

6)依据中各像素在图像graph[i]中的邻接关系，对中的像素进行分组排序，将在图像中依次邻接的像素按邻接顺序排列并归于一组，若只包含一组像素，执行步骤8)，否则执行步骤7)；

7)对中每个像素组的两个末端像素，分别找出与之距离最短的另一像素组的末端像素与之配对，以配对像素之一为起点，另一配对像素为终点，以起点和终点间的连线方向为检测方向，将逐次减1并检测空隙轮廓像素，直至起点和终点间的空隙轮廓像素在两者连线方向上依次相邻；

8)调用ACIS7.0的底层建模功能Spline_Fitting拟合中的像素，获得FACE类型的空隙轮廓面实体模型

9)使用空隙映射算法为各隶属于同一独立空隙的中的空隙编号赋值；

所述均匀方形网格划分是指从graph[i]中以一定间隔分别沿纵横方向抽取相应像素作为网格节点像素的网格划分方式；

所述已访问标志bV是指用以标识graph[i]中像素在空隙轮廓检测过程中访问状态的标志；

所述八邻域像素是指图像中一个像素的上、下、左、右、左上、右上、左下和右下八个方向上与之邻接的像素；

所述灰度差分八叉树是指用以描述graph[i]中第j个空隙轮廓在检测过程中所涉及像素间的邻接关系及邻接像素间的灰度差分的八叉树；

所述八叉树是指由若干结点按照结点间的父子关系以直线段连接而成、且每个结点最多有八个子结点的一种树形数据结构；

所述可扩展标志bE是指用以标识中的结点在的构建过程中能否作为父结点的标志；

所述灰度差分值d_FC是指赋予内子结点与父结点连接线段的、表示父子结点对应像素间灰度差的数值；

所述当前层号h是指在构建过程中的当前层数；

所述灰度差异阈值是指用以判断graph[i]中两任意相邻像素是否分别为空隙和沥青胶浆或集料像素中的任一类像素的灰度差分最小值，依据图像中局部区域内空隙像素间的灰度差异特点而预先设定；

所述空隙轮廓元是指包含graph[i]中第j个空隙轮廓的轮廓像素链表 及该轮廓隶属的空隙编号的数据结构；

所述轮廓像素链表是指用以存储graph[i]中第j个空隙轮廓的轮廓像素的链表型数据结构；

所述空隙编号是指graph[i]中第j个空隙轮廓所隶属单个空隙的编号；

所述预设最大层号H是指所建允许包含的最大层数；

所述空隙轮廓面实体模型是指以P_j中像素的拟合结果为边界的平面的实体模型。

所述的基于X-ray CT图像的沥青混合料空隙空间结构重构方法，其特征在于：所述空隙像素判别算法描述如下：

2)将集料像素灰度范围上限赋值为graph[i]对应灰度直方图中正方向上的最大灰度值；

3)令灰度变量使g₁自增，直至区间内的像素占试件全部像素的比例超过空隙体积百分比p_V；

4)令灰度变量使g₂自减，直至区间内的像素占试件全部像素的比例超过集料体积百分比p_G；

5)若待判像素的灰度g在内，执行步骤6)，否则待判像素不是空隙像素，算法结束；

6)由待判像素出发，以集料级配参数中的最大粒径为最大散射距离d_max进行四方向散射，任一方向的散射过程至散射距离等于d_max或检测到该散射方向上两个相邻像素间的灰度差异绝对值大于时结束，若所有方向上均进行了散射距离大于1的散射，散射结束的原因均为检测到该散射方向上两个相邻像素间的灰度差异是大于的正值，且两相邻像素中位于散射方向前端的像素的灰度不在 内，则待判像素为某一空隙轮廓内的空隙像素；

所述图像面积a_I是指一幅图像中全部像素的数量；

所述试件截面面积a_S是指图像中的沥青混合料试件像素的数量；

所述零灰度像素比例r_z是指图像中灰度为0的像素占图像全部像素的比例；

所述空隙像素灰度范围下限是指图像中空隙像素灰度范围内的灰度最小值；

所述灰度直方图是指相应图像中各灰度级上像素数量的统计图；

所述集料像素灰度范围上限是指图像中集料像素灰度范围内的灰度最大值；

所述空隙体积百分比p_V是指沥青混合料试件中空隙的体积比例；

所述集料体积百分比p_G是指沥青混合料试件中集料的体积比例；

所述最大散射距离d_max是指在四方向散射过程中任一方向上散射的最大距离；

所述四方向散射是指由起始像素出发，以一个像素为步长，在最大散射距离内分别沿上、下、左、右四个方向逐步散射，每散射一步，将散射终点像素的灰 度减去散射起点像素的灰度，用以判断散射是否结束的过程。

所述的基于X-ray CT图像的沥青混合料空隙空间结构重构方法，其特征在于：所述空隙映射算法描述如下：

1)令变量i＝1,j＝1；

3)判断中像素的灰度平均值与中像素的灰度平均值之间的相对误差是否小于灰度连续阈值T_C，若是，将和分别对应的和均赋值为j，否则仅将赋值为j；

4)使j自增1并执行步骤2)～步骤4)直至j＝q_i；

5)使i自增1并执行步骤2)～步骤5)直至i＝n；

所述包围盒是指将完全包围起来的最小长方形；

所述面积重合阈值T_A是指预设用于初步判定两个空隙轮廓面实体模型隶属于同一个空隙、相应两个包围盒的面积重合比例最小值；

所述灰度连续阈值T_C是指预设用于判定两个空隙轮廓面实体模型隶属于同一个空隙、相应两组轮廓像素灰度平均值的相对误差最大值。

所述的基于X-ray CT图像的沥青混合料空隙空间结构重构方法，其特征在于，所述的ACIS7.0是由美国Spatial公司生产的基于面向对象软件技术的三维几何造型引擎。

本发明的有益效果如下：

本发明基于沥青混合料细观结构中的空隙分布统计，完成从沥青混合料试件X-ray CT断面图像到空隙空间结构实体模型的自动转换，以辅助统计沥青混合 料内部空隙的真实空间分布，并为基于真实细观结构的沥青混合料细观三维数值建模奠定基础。

附图说明

图1为某沥青混合料试件X-ray CT扫描图像序列。

图2为序列各图像中空隙轮廓像素的检测结果图。

图3a为某单个空隙在序列各图像中对应的轮廓面实体模型。

图3b为某单个空隙对应所有轮廓面实体模型的边界。

图3c为某单个空隙对应的实体模型。

图4为沥青混合料试件中的空隙空间结构实体模型。

具体实施方式

基于X-ray CT图像的沥青混合料空隙空间结构重构方法，包括以下步骤：

(1)输入沥青混合料试件的X-ray CT断面扫描图像文件组graph[n]，用以重构空隙空间结构的实体模型；

沥青混合料试件是指以沥青材料和集料按适当比例配制，胶结成整体的圆柱形样本，记半径为radius，高度为height；

X-ray CT断面扫描图像文件组是指使用X-ray CT设备以一定纵向长度间隔扫描沥青混合料试件获取的图像数组，其中，图像数组中的元素均为BMP格式、分辨率为512×512像素的灰度图像文件；

灰度图像文件是指将彩色图像中像素的颜色分量值按照一定比例进行加权，并将加权值附于每个颜色分量而获得的图像文件，其中，灰度图像中所有像素的灰度均介于0到255之间；

空隙空间结构是指反映沥青混合料试件中所有空隙的几何形状及分布信息的空间结构；

实体模型是指以边界表示法表示的三维几何结构，记为模型M＝(FACE,EDGE,VERTEX,R)，其中FACE表示模型中面的集合，EDGE表示边的集合，VERTEX表示顶点的集合，R表示模型中各元素之间的邻接关系； 记录着面的类别：平面、样条面，以及面的曲面方程，其中face为FACE中的元素；

记录着边的曲线方程，其中edge为EDGE中的元素；

记录着顶点的坐标p(x,y,z)，其中vertex为VERTEX中的元素；

R＝{r₁,r₂}，其中：

r₁＝{(face₁,face₂,edge)|edge∈face₁∩face₂,edge∈EDGE,face_i∈FACE,i＝1,2}，r₁表示face₁和face₂相邻于边edge；

r₂＝{(edge₁,edge₂,vertex)|vertex∈edge₁∩edge₂,vertex∈VERTEX,edge_i∈EDGE,i＝1,2}，r₂表示edge₁和edge₂相邻于点vertex；

(2)使用空隙轮廓检测算法，识别出每幅X-ray CT断面扫描图像graph[i](i＝1,2,…,n)中各单个空隙的轮廓像素，获得相应的空隙轮廓面实体模型，并建立沥青混合料中各单个空隙与相应的轮廓面实体模型间的映射关系；

空隙轮廓检测算法描述如下：

1)对graph[i]进行纵横向间距均为3个像素的均匀方形网格划分，获得网格节点像素集合NS_i{node_p|p＝1,2,…,16384}，其中各网格节点像素按其在graph[i]中的位置自上而下、自左至右进行顺序编号，令p＝1、graph[i]中所有像素的已访问标志bV＝false；

2)使用空隙像素判别算法判断node_p是否为位于某一空隙轮廓内的空隙像素，若是，将node_p及其八邻域像素分别作为该空隙轮廓对应灰度差分八叉树 的根结点及根结点的子结点，置此9个像素的已访问标志bV＝true并置根结点的所有子结点的可扩展标志bE＝true后计算中每条连接线段的灰度差分值d_FC，令当前层号h＝2后执行步骤3)；若否，仅置node_p的已访问标志bV＝true，使p增1后执行步骤2)；

3)对第h层上的结点，依次判断每个结点与其父结点的连接线段对 应的d_FC是否超过灰度差异阈值若是，则该结点对应的像素为空隙轮廓像素，将该像素存入空隙轮廓元中的轮廓像素链表并设置该结点的可扩展标志bE＝false；

4)对第h层上bE＝true的结点，将其对应像素的八邻域像素中bV＝false的像素作为该结点的子结点并置各子结点所对应像素的bV＝true，使h增1；

5)循环执行步骤3)～步骤4)，直至h等于预设最大层号H或第h层上的结点所对应像素均为空隙轮廓像素，完成构建

6)依据中各像素在图像graph[i]中的邻接关系，对中的像素进行分组排序，将在图像中依次邻接的像素按邻接顺序排列并归于一组，若只包含一组像素，执行步骤8)，否则执行步骤7)；

7)对中每个像素组的两个末端像素，分别找出与之距离最短的另一像素组的末端像素与之配对，以配对像素之一为起点，另一配对像素为终点，以起点和终点间的连线方向为检测方向，将逐次减1并检测空隙轮廓像素，直至起点和终点间的空隙轮廓像素在两者连线方向上依次相邻；

8)调用ACIS7.0的底层建模功能Spline_Fitting拟合中的像素，获得FACE类型的空隙轮廓面实体模型

9)使用空隙映射算法为各隶属于同一独立空隙的中的空隙编号赋值；

均匀方形网格划分是指从graph[i]中以一定间隔分别沿纵横方向抽取相应像素作为网格节点像素的网格划分方式；

已访问标志bV是指用以标识graph[i]中像素在空隙轮廓检测过程中访问状态的标志；

空隙像素判别算法描述如下：

2)将集料像素灰度范围上限赋值为graph[i]对应灰度直方图中正方向上的最大灰度值；

3)令灰度变量使g₁自增，直至区间内的像素占试件全部像素的比例超过空隙体积百分比p_V；

4)令灰度变量使g₂自减，直至区间内的像素占试件全部像素的比例超过集料体积百分比p_G；

5)若待判像素的灰度g在内，执行步骤6)，否则待判像素不是空隙像素，算法结束；

6)由待判像素出发，以集料级配参数中的最大粒径为最大散射距离d_max进行四方向散射，任一方向的散射过程至散射距离等于d_max或检测到该散射方向上两个相邻像素间的灰度差异绝对值大于时结束，若所有方向上均进行了散射距离大于1的散射，散射结束的原因均为检测到该散射方向上两个相邻像素间的灰度差异是大于的正值，且两相邻像素中位于散射方向前端的像素的灰度不在 内，则待判像素为某一空隙轮廓内的空隙像素；

图像面积a_I是指一幅图像中全部像素的数量；

试件截面面积a_S是指图像中的沥青混合料试件像素的数量；

零灰度像素比例r_z是指图像中灰度为0的像素占图像全部像素的比例；

空隙像素灰度范围下限是指图像中空隙像素灰度范围内的灰度最小值；

灰度直方图是指相应图像中各灰度级上像素数量的统计图；

集料像素灰度范围上限是指图像中集料像素灰度范围内的灰度最大值；

空隙体积百分比p_V是指沥青混合料试件中空隙的体积比例；

集料体积百分比p_G是指沥青混合料试件中集料的体积比例；

最大散射距离d_max是指在四方向散射过程中任一方向上散射的最大距离；

四方向散射是指由起始像素出发，以一个像素为步长，在最大散射距离内分别沿上、下、左、右四个方向逐步散射，每散射一步，将散射终点像素的灰度减去散射起点像素的灰度，用以判断散射是否结束的过程；

八邻域像素是指图像中一个像素的上、下、左、右、左上、右上、左下和右下八个方向上与之邻接的像素；

灰度差分八叉树是指用以描述graph[i]中第j个空隙轮廓在检测过程中所涉及像素间的邻接关系及邻接像素间的灰度差分的八叉树；

八叉树是指由若干结点按照结点间的父子关系以直线段连接而成、且每个结点最多有八个子结点的一种树形数据结构；

可扩展标志bE是指用以标识中的结点在的构建过程中能否作为父结点的标志；

灰度差分值d_FC是指赋予内子结点与父结点连接线段的、表示父子结点对应像素间灰度差的数值；

当前层号h是指在构建过程中的当前层数；

灰度差异阈值是指用以判断graph[i]中两任意相邻像素是否分别为空隙和沥青胶浆或集料像素中的任一类像素的灰度差分最小值，依据图像中局部区域内空隙像素间的灰度差异特点而预先设定；

空隙轮廓元是指包含graph[i]中第j个空隙轮廓的轮廓像素链表及该轮廓隶属的空隙编号的数据结构；

轮廓像素链表是指用以存储graph[i]中第j个空隙轮廓的轮廓像素的链表型数据结构；

空隙编号是指graph[i]中第j个空隙轮廓所隶属单个空隙的编号；

预设最大层号H是指所建允许包含的最大层数；

空隙轮廓面实体模型是指以P_j中像素的拟合结果为边界的平面的实体模型；

空隙映射算法描述如下：

1)令变量i＝1,j＝1；

3)判断中像素的灰度平均值与中像素的灰度平均值之间的相对误差是否小于灰度连续阈值T_C，若是，将和分别对应的和均赋值为j，否则仅将赋值为j；

4)使j自增1并执行步骤2)～步骤4)直至j＝q_i；

5)使i自增1并执行步骤2)～步骤5)直至i＝n；

包围盒是指将完全包围起来的最小长方形；

面积重合阈值T_A是指预设用于初步判定两个空隙轮廓面实体模型隶属于同一个空隙、相应两个包围盒的面积重合比例最小值；

灰度连续阈值T_C是指预设用于判定两个空隙轮廓面实体模型隶属于同一个空隙、相应两组轮廓像素灰度平均值的相对误差最大值；

(3)输入各单个空隙所对应的全部轮廓面实体模型(k＝1,2,…,m)，调用ACIS7.0的底层建模功能Skinning，获得BODY类型的单个空隙实体模型M_k(k＝1,2,…,m)；

(4)输入各单个空隙实体模型{M_k|k＝1,2,…,m)}，调用ACIS7.0的底层建模功能Boolean_Unite，获得BODY类型的空隙空间结构实体模型M，该模型以边界表示法表示各单个空隙的三维几何结构和在沥青混合料试件中的分布规律，最终存储为SAT格式的图形文件，可应用于沥青混合料内部空隙的真实空间分布统计及沥青混合料细观数值分析；

空隙空间结构实体模型M是指能够反映沥青混合料试件中一定尺寸之上的所有空隙的形状属性及在试件中分布规律的空隙空间结构三维模型；

沥青混合料细观数值分析是指基于沥青混合料细观结构、用于预测沥青混合料的力学和热传导等性能的数值仿真。

本发明中，用C++语言，基于GDI+及ACIS内核，实现了本发明所描述的算法，并且以某沥青混合料试件的X-ray CT扫描图像为数据源，构建了试件内部空隙空间结构的实体模型。

(1)以直径为150mm，高度为164mm的圆柱形沥青混合料试件扫描步长为1mm的X-ray CT断面扫描图像文件组为数据源，输入如图1所示的X-ray CT断面扫描图像文件组graph[5]，用以重构空隙空间结构；

(2)为扫描图像文件graph[i](i＝1,2,3,4,5)设定灰度差异阈值如表1所示，并设定预设最大层号H＝15，空隙体积百分比p_V＝12.18％，集料体积百分比p_G＝67.84％，最大散射距离d_max＝15，面积重合阈值T_A＝0.8，灰度连续阈值T_C＝0.3，使用空隙轮廓检测算法，识别出graph[i](i＝1,2,3,4,5)中各单个空隙的轮廓像素结果如图2所示，获得相应的空隙轮廓面实体模型，并建立如图3a所示形式的沥青混合料中各单个空隙与相应的轮廓面实体模型间的映射关系；

表1各扫描图像文件对应的灰度差异阈值空隙轮廓元数量

(3)输入各单个空隙所对应的全部轮廓面实体模型(k＝1,2,…,m)，具有如图3b所示形式的空间排列顺序，调用ACIS7.0的底层建模功能Skinning，获得如图3c所示形式的BODY类型的单个空隙实体模型M_k(k＝1,2,…,m)；

(4)输入各单个空隙实体模型{M_k|k＝1,2,…,m)}，调用ACIS7.0的底层建模功能Boolean_Unite，获得如图4所示BODY类型的空隙空间结构实体模型M，该模型以边界表示法表示各单个空隙的三维几何结构和在沥青混合料试件中的分布规律，最终存储为SAT格式的图形文件，可应用于沥青混合料内部空隙的真实空间分布统计及沥青混合料细观数值分析。

Claims (5)

1.基于X-ray CT图像的沥青混合料空隙空间结构重构方法，其特征在于：包括以下步骤：

(1)输入沥青混合料试件的X-ray CT断面扫描图像文件组graph[n]，用以重构空隙空间结构的实体模型；

所述沥青混合料试件是指以沥青材料和集料按适当比例配制，胶结成整体的圆柱形样本，记半径为radius，高度为height；

所述X-ray CT断面扫描图像文件组是指使用X-ray CT设备以一定纵向长度间隔扫描沥青混合料试件获取的图像数组，其中，图像数组中的元素均为BMP格式，分辨率为512×512像素的灰度图像文件；

所述灰度图像文件是指将彩色图像中像素的颜色分量值按照一定比例进行加权，并将加权值附于每个颜色分量而获得的图像文件，其中，灰度图像中所有像素的灰度均介于0到255之间；

所述空隙空间结构是指反映沥青混合料试件中所有空隙的几何形状及分布信息的空间结构；

所述实体模型是指以边界表示法表示的三维几何结构，记为模型M＝(FACE,EDGE,VERTEX,R)，其中FACE表示模型中面的集合，EDGE表示边的集合，VERTEX表示顶点的集合，R表示模型中各元素之间的邻接关系；记录着面的类别：平面、样条面，以及面的曲面方程，其中face为FACE中的元素；

记录着边的曲线方程，其中edge为EDGE中的元素；

记录着顶点的坐标p(x,y,z)，其中vertex为VERTEX中的元素；

R＝{r₁,r₂}，其中：

r₁＝{(face₁,face₂,edge)|edge∈face₁∩face₂,edge∈EDGE,face_i∈FACE,i＝1,2}，r₁表示face₁和face₂相邻于边edge；

r₂＝{(edge₁,edge₂,vertex)|vertex∈edge₁∩edge₂,vertex∈VERTEX,edge_i∈EDGE,i＝1,2}，r₂表示edge₁和edge₂相邻于点vertex；

(2)使用空隙轮廓检测算法，识别出每幅X-ray CT断面扫描图像graph[i]中各单个空隙的轮廓像素，其中i＝1,2,…,n；获得相应的空隙轮廓面实体模型，并建立沥青混合料中各单个空隙与相应的轮廓面实体模型间的映射关系；

(3)输入各单个空隙所对应的全部轮廓面实体模型其中k＝1,2,…,m，调用ACIS7.0的底层建模功能Skinning，获得BODY类型的单个空隙实体模型M_k，其中k＝1,2,…，m；所述空隙轮廓面实体模型是指以P_j中像素的拟合结果为边界的平面的实体模型，所述轮廓像素链表是指用以存储graph[i]中第j个空隙轮廓的轮廓像素的链表型数据结构；

(4)输入各单个空隙实体模型M_k，其中k＝1,2,…，m；调用ACIS7.0的底层建模功能Boolean_Unite，获得BODY类型的空隙空间结构实体模型M，该模型以边界表示法表示各单个空隙的三维几何结构和在沥青混合料试件中的分布规律，最终存储为SAT格式的图形文件，可应用于沥青混合料内部空隙的真实空间分布统计及沥青混合料细观数值分析；

所述空隙空间结构实体模型M是指能够反映沥青混合料试件中一定尺寸之上的所有空隙的形状属性及在试件中分布规律的空隙空间结构三维模型；

所述沥青混合料细观数值分析是指基于沥青混合料细观结构，用于预测沥青混合料的力学和热传导等性能的数值仿真。

2.根据权利要求1所述的基于X-ray CT图像的沥青混合料空隙空间结构重构方法，其特征在于：步骤(2)中，空隙轮廓检测算法描述如下：

1)对graph[i]进行纵横向间距均为3个像素的均匀方形网格划分，获得网格节点像素集合NS_i{node_p|p＝1,2,…,16384}，其中各网格节点像素按其在graph[i]中的位置自上而下、自左至右进行顺序编号，令p＝1、graph[i]中所有像素的已访问标志bV＝false；

2)使用空隙像素判别算法判断node_p是否为位于某一空隙轮廓内的空隙像素，若是，将node_p及其八邻域像素分别作为该空隙轮廓对应灰度差分八叉树的根结点及根结点的子结点，置此9个像素的已访问标志bV＝true并置根结点的所有子结点的可扩展标志bE＝true后计算中每条连接线段的灰度差分值d_FC，令当前层号h＝2后执行步骤3)；若否，仅置node_p的已访问标志bV＝true，使p增1后执行步骤2)；

3)对第h层上的结点，依次判断每个结点与其父结点的连接线段对应的d_FC是否超过灰度差异阈值若是，则该结点对应的像素为空隙轮廓像素，将该像素存入空隙轮廓元中的轮廓像素链表并设置该结点的可扩展标志bE＝false；

4)对第h层上bE＝true的结点，将其对应像素的八邻域像素中bV＝false的像素作为该结点的子结点并置各子结点所对应像素的bV＝true，使h增1；

5)循环执行步骤3)～步骤4)，直至h等于预设最大层号H或第h层上的结点所对应像素均为空隙轮廓像素，完成构建

6)依据中各像素在图像graph[i]中的邻接关系，对中的像素进行分组排序，将在图像中依次邻接的像素按邻接顺序排列并归于一组，若只包含一组像素，执行步骤8)，否则执行步骤7)；

7)对中每个像素组的两个末端像素，分别找出与之距离最短的另一像素组的末端像素与之配对，以配对像素之一为起点，另一配对像素为终点，以起点和终点间的连线方向为检测方向，将逐次减1并检测空隙轮廓像素，直至起点和终点间的空隙轮廓像素在两者连线方向上依次相邻；

8)调用ACIS7.0的底层建模功能Spline_Fitting拟合中的像素，获得FACE类型的空隙轮廓面实体模型

9)使用空隙映射算法为各隶属于同一独立空隙的中的空隙编号赋值；

所述均匀方形网格划分是指从graph[i]中以一定间隔分别沿纵横方向抽取相应像素作为网格节点像素的网格划分方式；

所述已访问标志bV是指用以标识graph[i]中像素在空隙轮廓检测过程中访问状态的标志；

所述八邻域像素是指图像中一个像素的上、下、左、右、左上、右上、左下和右下八个方向上与之邻接的像素；

所述灰度差分八叉树是指用以描述graph[i]中第j个空隙轮廓在检测过程中所涉及像素间的邻接关系及邻接像素间的灰度差分的八叉树；

所述八叉树是指由若干结点按照结点间的父子关系以直线段连接而成、且每个结点最多有八个子结点的一种树形数据结构；

所述可扩展标志bE是指用以标识中的结点在的构建过程中能否作为父结点的标志；

所述灰度差分值d_FC是指赋予内子结点与父结点连接线段的、表示父子结点对应像素间灰度差的数值；

所述当前层号h是指在构建过程中的当前层数；

所述灰度差异阈值是指用以判断graph[i]中两任意相邻像素是否分别为空隙和沥青胶浆或集料像素中的任一类像素的灰度差分最小值，依据图像中局部区域内空隙像素间的灰度差异特点而预先设定；

所述空隙轮廓元是指包含graph[i]中第j个空隙轮廓的轮廓像素链表及该轮廓隶属的空隙编号的数据结构；

所述空隙编号是指graph[i]中第j个空隙轮廓所隶属单个空隙的编号；

所述预设最大层号H是指所建允许包含的最大层数。

3.根据权利要求2所述的基于X-ray CT图像的沥青混合料空隙空间结构重构方法，其特征在于：所述空隙像素判别算法描述如下：

1)将图像面积a_I减去试件截面面积a_S后除以a_I获得的数值与零灰度像素比例r_z相比较，若二者差值的绝对值小于3％，则将空隙像素灰度范围下限赋值为graph[i]对应灰度直方图中正方向上灰度大于0的首个像素的灰度值；否则令

2)将集料像素灰度范围上限赋值为graph[i]对应灰度直方图中正方向上的最大灰度值；

3)令灰度变量使g₁自增，直至区间内的像素占试件全部像素的比例超过空隙体积百分比p_V；

4)令灰度变量使g₂自减，直至区间内的像素占试件全部像素的比例超过集料体积百分比p_G；

5)若待判像素的灰度g在内，执行步骤6)，否则待判像素不是空隙像素，算法结束；

6)由待判像素出发，以集料级配参数中的最大粒径为最大散射距离d_max进行四方向散射，任一方向的散射过程至散射距离等于d_max或检测到该散射方向上两个相邻像素间的灰度差异绝对值大于时结束，若所有方向上均进行了散射距离大于1的散射，散射结束的原因均为检测到该散射方向上两个相邻像素间的灰度差异是大于的正值，且两相邻像素中位于散射方向前端的像素的灰度不在内，则待判像素为某一空隙轮廓内的空隙像素；

所述图像面积a_I是指一幅图像中全部像素的数量；

所述试件截面面积a_S是指图像中的沥青混合料试件像素的数量；

所述零灰度像素比例r_z是指图像中灰度为0的像素占图像全部像素的比例；

所述空隙像素灰度范围下限是指图像中空隙像素灰度范围内的灰度最小值；

所述灰度直方图是指相应图像中各灰度级上像素数量的统计图；

所述集料像素灰度范围上限是指图像中集料像素灰度范围内的灰度最大值；

所述空隙体积百分比p_V是指沥青混合料试件中空隙的体积比例；

所述集料体积百分比p_G是指沥青混合料试件中集料的体积比例；

所述最大散射距离d_max是指在四方向散射过程中任一方向上散射的最大距离；

所述四方向散射是指由起始像素出发，以一个像素为步长，在最大散射距离内分别沿上、下、左、右四个方向逐步散射，每散射一步，将散射终点像素的灰度减去散射起点像素的灰度，用以判断散射是否结束的过程。

4.根据权利要求2所述的基于X-ray CT图像的沥青混合料空隙空间结构重构方法，其特征在于：所述空隙映射算法描述如下：

1)令变量i＝1,j＝1；

2)取得的包围盒及graph[i]、graph[i+1]分别对应的空隙轮廓面实体模型数量q_i和q_i+1，在从graph[i+1]中检测到的空隙轮廓面实体模型集合{OF_s ⁽ⁱ⁺¹⁾|s＝1,2,…,q_i+1}中寻找满足大于面积重合阈值T_A的OF_s ⁽ⁱ⁺¹⁾，若存在相应OF_s ⁽ⁱ⁺¹⁾，执行步骤3)，若无相应OF_s ⁽ⁱ⁺¹⁾，令对应的使j自增1后执行步骤2)；

3)判断中像素的灰度平均值与P_s ⁽ⁱ⁺¹⁾中像素的灰度平均值之间的相对误差是否小于灰度连续阈值T_C，若是，将和OF_s ⁽ⁱ⁺¹⁾分别对应的和均赋值为j，否则仅将赋值为j；

4)使j自增1并执行步骤2)～步骤4)直至j＝q_i；

5)使i自增1并执行步骤2)～步骤5)直至i＝n；

所述包围盒是指将完全包围起来的最小长方形；

所述面积重合阈值T_A是指预设用于初步判定两个空隙轮廓面实体模型隶属于同一个空隙、相应两个包围盒的面积重合比例最小值；

所述灰度连续阈值T_C是指预设用于判定两个空隙轮廓面实体模型隶属于同一个空隙、相应两组轮廓像素灰度平均值的相对误差最大值。

5.根据权利要求1所述的基于X-ray CT图像的沥青混合料空隙空间结构重构方法，其特征在于，所述的ACIS7.0是由美国Spatial公司生产的基于面向对象软件技术的三维几何造型引擎。