一种航天器多类型复杂缺陷的量化评估方法与流程

1.本发明属于航空航天飞行器损伤检测与维修保障技术领域，更具体地说，本发明涉及一种航天器多类型复杂缺陷的量化评估方法


背景技术：

2.针对空间碎片与航天器发生超高速撞击形成的多类型缺陷的检测问题，目前已开展了很多基于红外热成像检测技术的缺陷量化评估方面的研究。利用红外热成像检测技术对缺陷进行检测，对红外热图像序列处理获得表征缺陷特征的红外重构图像后，能够通过定性分析初步确定缺陷的形态特征。为了保障航天器运行安全，需要进一步对缺陷特征进行的量化评估。基于航天核心构件的复杂多类型损伤缺陷特征的量化评估结果能够指导进行超高速撞击损伤研究、航天器防护结构设计、航天器损伤修护等，为航天器在轨安全服役提供保障。
3.超高速撞击损伤具有种类多样、类型复杂和散布区域广的特点给缺陷检测与量化评估带来了困难。并且在检测过程中，航天器表面材料飞溅的颗粒物和噪声干扰均会影响缺陷的定量检测结果。对于超高速撞击损伤检测困难和客观存在的检测环境影响问题，为了更好的对复杂多类型的超高速撞击损伤进行识别，将定量检测分为定位和量化两个部分进行。对不同类型和分布在不同位置的缺陷分别定位，对复杂多类型的缺陷特征识别和修复具有重要意义。根据定位结果分别量化不同类型和不同区域的缺陷特征信息，半年内、为试件受损情况提供具体的数值参考。从定位和量化两个方面来描述超高速撞击形成的复杂多类型缺陷的形态特征和物理特征，根据不同损伤的特征制定维修方案。
4.在2021年6月1日公布的、公布号为cn112881467a、名称为“一种大尺寸复合材料损伤成像与定量识别方法”的中国发明专利申请中，在对大尺寸复合材料进行损伤定量检测的过程中，通过分割的方法提取红外重构图像中的缺陷特征，并利用外接矩形来标定缺陷区域，再根据矩形框的标定结果分别对各个缺陷区域进行检测并量化缺陷形态特征。
5.通过外接矩形框标定缺陷区域的方法针对大尺寸复合材料能够很好的识别检测出缺陷。然而在应用于超高速撞击损伤定量检测的过程中，超高速撞击损伤类型复杂，通过外接矩形来标定的缺陷区域与实际的缺陷区域位置可能存在较大的偏差，只能对缺陷进行粗略的定位，会影响缺陷特征识别结果和缺陷修复的效率。
6.针对外接矩形框只能进行粗略定位的问题，可以通过边界追踪法提取缺陷边界轮廓。缺陷的边界轮廓相较于外接矩形框包含更多的缺陷细节信息，基于缺陷的边界轮廓来定位缺陷能够提高缺陷定位精度。但是像素之间存在不可细分的间隙，当缺陷边界落在像素之间的间隙处时，整像素的定位方法会出现较大的误差。并且离散的缺陷像素轮廓放大后会失真，影响缺陷区域判定。
7.并且在量化缺陷物理特征时，基于缺陷图像计算得到的结果实际是像素尺度上的数值特征，需要根据图像尺寸和待测试件实际尺寸的比例关系进行换算才能得到实际尺寸。然而在实际计算过程中，由于图像大小和待测试件实际尺寸不同，位于水平方向和垂直
方向的像素点的像素尺度不同，因此没有统一的比例来计算缺陷的实际物理信息，导致一些衡量缺陷物理特征的参数无法计算。
8.针对上述问题，本发明提出了一种航天器多类型复杂缺陷的量化评估方法。由于红外重构图像中的像素点反映连续的温度变化，没有明确的类别划分，对此本发明选择了模糊c均值算法来进行分割处理。考虑到缺陷类型复杂和客观存在的热扩散现象对初步提取的缺陷特征的影响，在定量检测前，对初步提取的缺陷特征图像进行了后处理，并分别对不同类型和不同位置的缺陷进行了区域修正。对处理后的缺陷特征图像，先分别提取不同类型、不同位置损伤的边界特征，对缺陷区域进行标定，实现缺陷像素粗定位。再对每一个标定的缺陷区域利用曲线拟合的方式实现亚像素缺陷特征定位。在量化各个缺陷区域的特征参数前，将拟合得到的缺陷轮廓进行尺寸映射，再进一步量化计算缺陷的面积、周长、长径和短径等特征。结合空间位置信息和物理特征量化信息来描述缺陷特征，以满足对多类型复杂缺陷进行定量识别的需求。


技术实现要素：

9.本发明的一个目的是解决至少上述问题和/或缺陷，并提供至少后面将说明的优点。
10.为了实现根据本发明的这些目的和其它优点，提供了一种航天器多类型复杂缺陷的量化评估方法，包括以下步骤：
11.步骤一、采集红外数据，对红外热图像序列数据处理，获得表征缺陷特征的红外重构图像igs，并利用模糊c聚类分割的方式对表征缺陷特征的红外重构图像igs进行缺陷特征提取得到缺陷特征提取图像ige；
12.步骤二、利用基于隶属度对缺陷特征提取图像ige进行后处理，优化分割结果，得到缺陷特征检测图像igt；
13.步骤三、通过时序排齐法对缺陷特征检测图像igt中的不同类型和不同位置的缺陷区域进行修正，得到缺陷特征识别图像igc；
14.步骤四、通过8-邻接freeman链码来检测缺陷特征识别图像igc中不同类型和不同位置缺陷的边界特征，对缺陷区域进行标定，实现缺陷像素粗定位；
15.步骤五、基于缺陷像素粗定位结果，根据提取到的缺陷边界像素点，对每一个标定的缺陷区域利用b
é
zier曲线拟合的方式实现亚像素缺陷特征定位；
16.步骤六、根据亚像素缺陷定位结果，对每一个标定的缺陷区域，对轮廓特征尺寸在水平方向和垂直方向，按照试件和图像的尺寸关系，分别进行映射以得到刻画缺陷不同特征参数的实际尺寸；
17.步骤七、进行缺陷特征量化计算。
18.优选的是，其中，所述步骤一中，进行缺陷特征提取得到缺陷特征提取图像ige的具体步骤包括：
19.步骤1.1、设置聚类数目为c，初始化聚类中心vk，k＝1,
…
,c，初始化缺陷滤波图像ige中像素点xa与聚类中心的隶属度vk的隶属度u
ak
，像素点xa的属度u
ak
满足
20.步骤1.2、计算并更新聚类中心vk：
[0021][0022]
其中，m为平滑参数，n为缺陷滤波图像igs中的像素点数目；
[0023]
步骤1.3、更新像素点xa的隶属度函数u
ak
：
[0024][0025]
步骤1.4、计算当前第t次迭代时的目标函数j
t
：
[0026][0027]
步骤1.5、迭代终止条件判断：如果第t次迭代时的目标函数j
t
与第t-1次迭代时的目标函数j
t-1
的差值小于迭代条件判断阈值ε，即||j
t-j
t-1
||＜ε，或者达到最大迭代次数t
max
，迭代终止，得到最终的聚类中心集合vk，k＝1,
…
,c和缺陷滤波图像igs各像素点的隶属度矩阵u
ak
；
[0028]
步骤1.6、缺陷滤波图像igs像素点xa划分为隶属度矩阵u
ak
值最大的那一类，缺陷滤波图像igs中所有的像素点xa划分结束得到最终的缺陷特征提取图像ige。
[0029]
优选的是，其中，所述步骤二的具体步骤包括：
[0030]
步骤2.1、检测所有可能的错分像素点；检测所有像素点xa的r
×
r邻域区域中，各邻域像素点的分类标签，如果分类标签相同，则保持原分类标签不变，反之则认为标签不一致的像素点可能出现错分的情况，将其提取保存；
[0031]
步骤2.2、对提取的可能错分的像素点进行标签重置；对可能错分的像素点xe，统计xe的n
×
n邻域区域中的像素点x
p
的隶属度，重新计算像素点xe属于第k类的隶属度u
ek
'：
[0032][0033]
其中，ωe表示像素点xe的r
×
r邻域像素点集合；表示xe的邻域像素点属于第k类的隶属度；γe表示ωe中像素点的数目，根据u
ek
'的计算情况，将xe划分为隶属度值最大的那一类，实现像素点标签重置；在进行像素点标签重置时需要注意，重置的标签必须是该像素点在步骤2.2中r
×
r邻域区域出现过的分类标签；对缺陷特征提取图像ige中所有可能出现的缺陷区域像素点进行分类标签重置后得到缺陷特征检测图像igt。
[0034]
优选的是，其中，所述步骤三中，将步骤二得到的缺陷特征检测图像igt分为k个区域，对第k，k＝1,
…
,k个缺陷区域，通过对缺陷区域像素点进行相似性度量判断以获取热扩散区域的大小：
[0035]
kdwei
≥kε
ttr_i
×
kdwemax
[0036]
其中，
kdwei
，表示第k个缺陷区域中，像素点kxi与缺陷特征实际区域参
考点之间的距离；为第k个缺陷区域的像素点数目；
kdwemax
为第k个缺陷区域中的像素点与缺陷特征实际区域参考点之间的最大距离；kε
ttr_i
为缺陷像素点kxi的热扩散判定阈值；通过分别度量缺陷区域像素点与缺陷特征实际区域参考点之间的温度热响应相关性和缺陷区域像素点与背景区域参考点之间的温度热响应相关性来计算第k个缺陷区域像素点kxi的热扩散判断阈值kε
ttr_i
；考虑到客观存在的热扩散现象对缺陷区域进行修正，得到缺陷特征识别图像igc的具体步骤为：
[0037]
步骤3.1、时序排齐法计算第k个缺陷区域像素点kxi的热扩散判断阈值kε
ttr_i
[0038]
步骤3.1.1、计算缺陷区域像素点kxi和背景区域参考点之间的瞬态热响应序列相关性系数；
[0039]
步骤3.1.1.1、判断需要进行时序移动的瞬态热响应序列；选择背景区域的特征瞬态热响应对应的像素点作为背景区域参考点其对应的温度热响应序列为bg
re
(bg
re1
,
…
,bg
ret
)，t为热序列采集时间；缺陷区域像素点kxi的瞬态热响应序列分别找到背景区域参考点的温度峰值所在帧数和缺陷区域像素点kxi的温度峰值所在帧数kt1；如果则是需要进行移动的瞬态热响应序列，记为ttr2，记bg
re
(bg
re1
,
…
,bg
ret
)为ttr1；否则，bg
re
(bg
re1
,
…
,bg
ret
)是需要进行移动的瞬态热响应序列，记为ttr2，记
[0040]
步骤3.1.1.2、记ttr1的温度峰值所在帧数为t
max
，设置每次移动tt个时间单位，循环参数tm设置为1，有对ttr2进行时序移动，每次移动tt个时间单位，将ttr2的子序列和ttr1的子序列对齐并计算两个它们的相关性系数
[0041][0042]
其中，e(
·
)表示瞬态热响应序列的期望；
[0043]
步骤3.1.1.3、tm＝tm+1，循环步骤3.1.1.2，每移动一个时间间隔，计算一次相应的皮尔逊相关系数，直到为止；得到相关系数集合
[0044]
步骤3.1.1.4、从中找到最大值作为第k个缺陷区域像素点kxi与背景区域参考点温度热响应序列相关系数kρ
i&bg
：
[0045][0046]
步骤3.1.2、按照步骤3.1.1相同的计算方式，计算第k个缺陷区域像素点与缺陷特
征实际区域参考点温度热响应的相关系数kρ
i&def
；根据相关系数kρ
i&bg
和kρ
i&def
得到第k个缺陷区域像素点xi的热扩散区域判定阈值kε
ttr_i
：
[0047][0048]
步骤3.2、根据计算得到的热扩散区域判定阈值kε
ttr_i
，当缺陷像素点与缺陷特征实际区域参考点之间的距离
kdwei
满足：
[0049][0050]
认为像素点xi位于热扩散区域，则去掉像素点xi；重复步骤3.2直到缺陷特征识别图像igc中第k个缺陷区域像素点xi都被访问到，去除所有位于热扩散区域的像素点xi，转至步骤3.3；
[0051]
步骤3.3、重复步骤3.1和步骤3.2直到k个缺陷区域像素点均被访问到，得到缺陷特征识别图像igc。
[0052]
优选的是，其中，所述步骤四通过8-邻接freeman链码来检测缺陷特征识别图像igc中不同类型和不同位置缺陷的边界特征，对缺陷区域进行标定，实现缺陷像素粗定位，其中，对第k，k＝1,
…
,k个缺陷区域进行像素标定定位的具体过程为：
[0053]
选择第k个缺陷区域的边界跟踪的起始像素点，记为ke0，以链码值dir＝0的方向为开始边界搜索的方向，逆时针方向旋转45
°
，按照方向dir＝1搜索下一个像素点；由于经过提取得到的缺陷特征图像为二值图像，在这个二值图像中，表征缺陷轮廓像素点的像素值用“1”表示；基于此，在缺陷边界像素点的搜索过程中，如果像素点kx
t
的像素值为“1”，则认为像素点kx
t
是缺陷边界点，将其记为
ket
并保存；并将kx
t
的链码值赋值给它的前一个点像素kx
t-1
；将像素点kx
t
作为新的边界搜索起始点，逆时针方向旋转90
°
，按照方向dir＝2搜索下一个缺陷边界点；重复上述步骤，当搜索到起始点kb0时，缺陷区域搜索完毕；整个搜索过程中，检测到的边界点ke0,
…
,
kem
即为缺陷的边界轮廓；边界点集合kω
l
，kω
l
＝(ke0,
…
,
kem
)对缺陷进行了像素粗定位；重复上述过程直到k个缺陷区域访问完毕，得到k个缺陷区域的像素标定定位结果。
[0054]
优选的是，其中，所述步骤五的具体步骤包括：通过像素粗定位后，标定得到了k个缺陷区域，提取标定的第k，k＝1,
…
,k个缺陷区域的m个缺陷边界像素点，m个缺陷边界像素点构成像素点集合
kei
，i＝0,1,
…
,m；分别对相邻的缺陷边界像素点
kei
和像素点
kei+1
两两进行三次b
é
zier曲线拟合；由于i的最大值为m，因此像素点
kem+1
是不存在的；考虑到缺陷轮廓是封闭图形，用ke0作为
kem+1
与
kem
进行三次b
é
zier曲线拟合；利用三次b
é
zier曲线对缺陷边界像素点
kei
和
kei+1
两两进行拟合得到亚像素缺陷封闭轮廓的具体过程为：
[0055]
步骤5.1、在待拟合的边缘轮廓像素点
kei
和
kei+1
之间，构造控制点kki和
kci
：
[0056]
[0057][0058]
其中，
ka*
和
kb*
是任意给定的正数；(kxi,kyi)是缺陷边界像素点
kei
的位置参数；由于经过边界追踪算法检测得到缺陷的m个边界像素点
kei
(kxi,kyi)，i＝0,1,
…
,m，因此缺陷边界像素点集合中不存在像素点ke-1
、
kem+1
、
kem+2
，其对应的位置参数(kx-1
,ky-1
)，(kx
m+1
,ky
m+1
)，(kx
m+2
,ky
m+2
)也就不存在；考虑到缺陷轮廓是闭合的，用(kxm,kym)的值作为(kx-1
,ky-1
)的值，(kx0,ky0)的值作为(kx
m+1
,ky
m+1
)的值，(kx1,ky1)的值作为(kx
m+2
,ky
m+2
)的值；
[0059]
步骤5.2、根据b
é
zier曲线的定义，缺陷边界像素点
kei
和
kei+1
之间的三次b
é
zier曲线的表达式为：
[0060][0061]
根据三次b
é
zier曲线的表达式得到缺陷边界像素点
kei
和像素点
kei+1
之间的拟合曲线；
[0062]
步骤5.3、步骤四遍历缺陷边界像素点集合
kei
，i＝0,1,
…
,m所有的像素点，分别进行两两拟合得到相邻像素点之间的拟合曲线；当i＝m时，用ke0作为
kem+1
与
kem
进行三次b
é
zier曲线拟合；
[0063]
步骤5.4、所有缺陷边界像素点遍历结束，得到m条三次b
é
zier曲线；m条b
é
zier曲线平滑连接得到第k个缺陷区域的亚像素缺陷封闭轮廓kl；基于封闭轮廓kl实现缺陷亚像素定位；
[0064]
步骤5.5、重复上述过程直到k个缺陷区域均被访问到，得到不同类型和不同位置的缺陷亚像素定位结果。
[0065]
优选的是，其中，所述步骤六映射得到刻划缺陷不同特征参数的实际尺寸的方法包括：对标定的第k，k＝1,
…
,k个缺陷区域，亚像素边界点的坐标记为(kxs,kys)，按照如下公式进行映射：
[0066][0067]
式中，l
×
w为检测试件的实际尺寸；m
×
n为拍摄得到的红外图像的尺寸，得到映射后的亚像素边界点坐标(kxs',kys')。
[0068]
优选的是，其中，所述步骤七进行缺陷特征量化计算的具体步骤包括：周长和面积是常规且直接地反映缺陷的尺寸信息，基于几何结构特征参数数值，可直接量化评估各个缺陷特征区域的尺寸大小信息，进而体现各个缺陷特征区域之间的几何尺寸差异；缺陷的面积可以用来量化缺陷在图像中的大小；分别计算不同的缺陷区域的特征参数；对标定的第k，k＝1,
…
,k个缺陷区域，对将映射后的亚像素边界点(kxs',kys')集合勾画的缺陷实际尺寸中的边界轮廓记为kl；由于大多数缺陷特征区域呈现不规则性，因此缺陷边界轮廓kl的形状呈不规则图形；直接计算缺陷边界轮廓kl包围的缺陷区域的面积比较困难；对此，用缺陷边界轮廓kl包围的区域kω的像素点kpi的数目和作为面积特征参数；由于在步骤七中进行了边界尺寸映射，因此不再需要进行试件与图像之间的面积换算，此时统计得到缺陷边界轮廓kl包围的区域kω的像素点数目即为缺陷的面积特征参数：
[0069][0070]
通过统计缺陷轮廓包围的区域面积能够进一步确定缺陷的位置信息，达到标定缺陷的目的；在经映射后的缺陷轮廓图像中，基于缺陷边界轮廓特征进一步计算缺陷的周长、长径和短径特征；在计算缺陷的周长时，由于拟合得到的缺陷亚像素轮廓是通过按顺序两两拟合缺陷边界像素点
kei
和
kei+1
得到的，并且缺陷边界特征尺寸映射不会改变亚像素边界点之间的位置关系，因此通过曲线积分可以求得边界像素点
kei
和
kei+1
之间的拟合曲线的实际长度kpi：
[0071][0072]
其中，为缺陷边界像素点
kei
和
kei+1
之间的三次b
é
zier曲线的表达式；m个边界像素点两两之间的实际拟合曲线长度的和即为缺陷的周长kp：
[0073][0074]
同时由于缺陷特征区域的不规则性，基于亚像素缺陷边界轮廓得到的最小外接矩形框的尺寸同样也可以反映缺陷的形态特征；根据亚像素缺陷边界轮廓的拟合轮廓kl，找到其最小外接矩形，用最小外接矩形来刻画缺陷的分布情况；最小外接矩形的长边即为缺陷的长径kdl，短边即为缺陷的短径kds；通过面积来利用等间隔旋转法来求解亚像素缺陷边界轮廓的最小外接矩形：亚像素缺陷边界轮廓图像在90
°
范围内等间隔地旋转，每次旋转都需要记录亚像素缺陷边界轮廓的外接矩形参数：分别找到亚像素缺陷边界轮廓在水平方向的坐标最大值kx
max
和最小值kx
min
，竖直方向的坐标最大值ky
max
和最小值ky
min
；则标定的第k个缺陷区域的外接矩形的长边为：
[0075]kdl＝max(ky
max-kymin
,kx
max-k
x
min
)
[0076]
短边为：
[0077]kds＝min(ky
max-kymin
,kx
max-k
x
min
)；
[0078]
则外接矩形的面积为kdl
×kds，本发明通过面积来定义最小外接矩形，计算最小外接矩形的具体过程为：首先根据亚像素缺陷边界轮廓分布坐标的最大值最小值能够得到一个没有旋转角度的亚像素缺陷边界轮廓区域的外接矩形，并计算当前状态下的最小外接矩形的长边、短边和面积，此时旋转角度α＝0
°
；设置等间隔旋转角度θ，对亚像素缺陷边界轮廓进行旋转一个角度θ，α＝α+θ，按照相同的方式求取亚像素缺陷边界轮廓的外接矩形，并记录相关参数；重复旋转过程，直到旋转角度α＝90
°
。最后，根据面积大小确定最小外接矩形；根据上述量化过程，基于亚像素缺陷边界轮廓拟合结果计算得到缺陷的物理特征参数。
[0079]
本发明至少包括以下有益效果：
[0080]
(1)、本发明提出了一种航天器多类型复杂缺陷的量化评估方法，将定量检测分为定位和量化两个环节，结合空间位置信息和物理特征量化信息来描述缺陷特征。在定位过程中，考虑到缺陷类型复杂和噪声的存在会干扰初步特征提取结果，利用了基于隶属度的模糊聚类结果重新分类的方法，对可能出现错分的像素点重新进行划分。在利用红外检测技术进行损伤检测时，由于热扩散现象是客观存在的，为了消除热扩散对缺陷识别的影响，
分别计算了缺陷区域像素点瞬态热响应与缺陷区域参考瞬态热响应的相关系数和背景区域参考瞬态热响应的相关系数，根据相关系数来自动确定扩散区域判定阈值，利用判定阈值来去除热扩散区域。根据区域修正结果，通过边界追踪的方法提取缺陷边界特征来对缺陷像素进行定位。然后，基于三次b
é
zier曲线亚像素轮廓拟合的方法来对缺陷准确定位。在量化部分，通过缺陷边界特征尺寸映射解决像素尺度不同的问题，计算了缺陷的面积、周长、长径和短径等参数，从多个方面来描述缺陷特征。
[0081]
(2)、本发明通过隶属度对缺陷特征提取图像进行后处理和时序排齐的方法修正缺陷区域，能够对基于分割方式的缺陷特征提取结果进行优化，并且能够去除红外检测技术中客观存在的热扩散现象的影响，有效提高缺陷定位和量化的准确率。
[0082]
(3)、本发明通过亚像素拟合缺陷轮廓的方式定位缺陷区域，亚像素能够将像素之间微观上存在的间隙进行细分，提高缺陷定位的精度。准确的缺陷定位结果不但能够帮助快速寻找缺陷形态和位置，还能够指导进行微小缺陷的二次局部检测，从而提高缺陷的检测效率。并且拟合得到的缺陷轮廓连续光滑，放大后像素块会出现失真的情况，而拟合曲线仍然保持连续，具有稳定的缺陷特征。
[0083]
(4)、本发明在量化缺陷特征前进行了缺陷边界特征尺寸映射，解决了像素尺度不同的问题。映射后，像素在水平方向和竖直方向的像素尺度和实际尺度是相同的，通过定量计算得到的参数尺寸即为缺陷特征的实际尺寸，结合不同的几何特征信息来对缺陷的全局特征进行表征，为试件受损情况提供具体的数值参考。
[0084]
本发明的其它优点、目标和特征将部分通过下面的说明体现，部分还将通过对本发明的研究和实践而为本领域的技术人员所理解。
附图说明
[0085]
图1为本发明提供的航天器多类型复杂缺陷的量化评估方法流程图；
[0086]
图2为本发明实施例中表征缺陷特征的红外重构图像；
[0087]
图3为本发明实施例中缺陷特征提取图像；
[0088]
图4为本发明实施例中缺陷特征检测图像；
[0089]
图5为本发明实施例中缺陷区域标记图像；
[0090]
图6为本发明实施例中缺陷特征识别图像；
[0091]
图7为本发明实施例中缺陷区域边界提取图像；
[0092]
图8为本发明实施例中亚像素缺陷区域定位图像。
具体实施方式
[0093]
下面结合附图对本发明做进一步的详细说明，以令本领域技术人员参照说明书文字能够据以实施。
[0094]
应当理解，本文所使用的诸如“具有”、“包含”以及“包括”术语并不排除一个或多个其它元件或其组合的存在或添加。
[0095]
本发明的目的在于克服现有技术的不足，提供了一种航天器多类型复杂缺陷的量化评估方法。在定量计算缺陷特征前，采用分割的方式提取缺陷特征，并考虑到复杂多类型的缺陷会影响特征提取的准确率和客观存在的热扩散效应会干扰定量检测结果，对分割初
步提取的缺陷特征图像进行了后处理，再分别对不同类型和不同位置的缺陷区域进行修正。在前期对大尺寸复合材料损伤进行定量检测的基础上，为了应对航天器超高速撞击形成的多类型复杂缺陷的定量检测需求，先通过提取缺陷边界来标定不同类型和不同损伤的缺陷特征，进一步再基于像素定位结果，对每一个标定的缺陷区域利用曲线拟合的方式实现亚像素缺陷特征定位。分别对拟合得到的亚像素轮廓进行尺寸映射，再量化缺陷特征尺寸参数。结合空间位置信息和物理特征量化信息来描述缺陷特征，以满足对多类型复杂缺陷进行定量识别的需求。
[0096]
如图1所示，为实现上述发明目的，本发明提出一种航天器多类型复杂缺陷的量化评估方法，其特征在于，包括以下步骤：
[0097]
步骤一、采集红外数据，对红外热图像序列数据处理，获得表征缺陷特征的红外重构图像igs。并利用模糊c聚类分割的方式对表征缺陷特征的红外重构图像igs进行缺陷特征提取得到缺陷特征提取图像ige，解决红外重构图像反映的是温度变化，各区域之间的划分并不明确，不能为视觉系统所直接使用的问题。缺陷特征提取的具体步骤为：
[0098]
步骤1.1、设置聚类数目为c，初始化聚类中心vk，k＝1,
…
,c，初始化缺陷滤波图像ige中像素点xa与聚类中心的隶属度vk的隶属度u
ak
，像素点xa的属度u
ak
满足
[0099]
步骤1.2、计算并更新聚类中心vk：
[0100][0101]
其中，m为平滑参数，n为缺陷滤波图像igs中的像素点数目。
[0102]
步骤1.3、更新像素点xa的隶属度函数u
ak
：
[0103][0104]
步骤1.4、计算当前第t次迭代时的目标函数j
t
[0105][0106]
步骤1.5、迭代终止条件判断：如果第t次迭代时的目标函数j
t
与第t-1次迭代时的目标函数j
t-1
的差值小于迭代条件判断阈值ε，即||j
t-j
t-1
||＜ε，或者达到最大迭代次数t
max
，迭代终止。得到最终的聚类中心集合vk，k＝1,
…
,c和缺陷滤波图像igs各像素点的隶属度矩阵u
ak
[0107]
步骤1.6、缺陷滤波图像igs像素点xa划分为隶属度矩阵u
ak
值最大的那一类。缺陷滤波图像igs中所有的像素点xa划分结束得到最终的缺陷特征提取图像ige。
[0108]
步骤二、利用基于隶属度对缺陷特征提取图像ige进行后处理，优化分割结果，得到缺陷特征检测图像igt。具体步骤为：
[0109]
步骤2.1、检测所有可能的错分像素点。检测所有像素点xa的r
×
r邻域区域中，各
邻域像素点的分类标签，如果分类标签相同，则保持原分类标签不变，反之则认为标签不一致的像素点可能出现错分的情况，将其提取保存。
[0110]
步骤2.2、对提取的可能错分的像素点进行标签重置。对可能错分的像素点xe，统计xe的n
×
n邻域区域中的像素点x
p
的隶属度，重新计算像素点xe属于第k类的隶属度u
ek
'：
[0111][0112]
其中，ωe表示像素点xe的r
×
r邻域像素点集合；表示xe的邻域像素点属于第k类的隶属度；γe表示ωe中像素点的数目。根据u
ek
'的计算情况，将xe划分为隶属度值最大的那一类，实现像素点标签重置。在进行像素点标签重置时需要注意，重置的标签必须是该像素点在步骤(2.2)中r
×
r邻域区域出现过的分类标签。对缺陷特征提取图像ige中所有可能出现的缺陷区域像素点进行分类标签重置后得到缺陷特征检测图像igt。
[0113]
步骤三、通过时序排齐法对缺陷特征检测图像igt中的不同类型和不同位置的缺陷区域进行修正，得到缺陷特征识别图像igc。经特征提取后，将缺陷分为了k个区域，对第k，k＝1,
…
,k个缺陷区域，通过对缺陷区域像素点进行相似性度量判断以获取热扩散区域的大小：
[0114]
kdwei
≥kε
ttr_i
×
kdwemax
[0115]
其中，
kdwei
，表示第k个缺陷区域中，像素点kxi与缺陷特征实际区域参考点之间的距离；为第k个缺陷区域的像素点数目；
kdwemax
为第k个缺陷区域中的像素点与缺陷特征实际区域参考点之间的最大距离；kε
ttr_i
为缺陷像素点kxi的热扩散判定阈值。通过分别度量缺陷区域像素点与缺陷特征实际区域参考点之间的温度热响应相关性和缺陷区域像素点与背景区域参考点之间的温度热响应相关性来计算第k个缺陷区域像素点kxi的热扩散判断阈值kε
ttr_i
。考虑到客观存在的热扩散现象对缺陷区域进行修正，得到缺陷特征识别图像igc的具体步骤为：
[0116]
步骤3.1、时序排齐法计算第k个缺陷区域像素点kxi的热扩散判断阈值kε
ttr_i
[0117]
步骤3.1.1、计算缺陷区域像素点kxi和背景区域参考点之间的瞬态热响应序列相关性系数。
[0118]
步骤3.1.1.1、判断需要进行时序移动的瞬态热响应序列。选择背景区域的特征瞬态热响应对应的像素点作为背景区域参考点其对应的温度热响应序列为bg
re
(bg
re1
,
…
,bg
ret
)，t为热序列采集时间。缺陷区域像素点kxi的瞬态热响应序列分别找到背景区域参考点的温度峰值所在帧数和缺陷区域像素点kxi的温度峰值所在帧数kt1。如果则是需要进行移动的瞬态热响应序列，记为ttr2，记bg
re
(bg
re1
,
…
,bg
ret
)为ttr1。否则，bg
re
(bg
re1
,
…
,bg
ret
)是需要进行移动的瞬态热响应序列，记为ttr2，记为ttr1。
[0119]
步骤3.1.1.2、记ttr1的温度峰值所在帧数为t
max
，设置每次移动tt个时间单位，循
环参数tm设置为1，有对ttr2进行时序移动，每次移动tt个时间单位，将ttr2的子序列和ttr1的子序列对齐并计算两个它们的相关性系数
[0120][0121]
其中，e(
·
)表示瞬态热响应序列的期望。
[0122]
步骤3.1.1.3、tm＝tm+1，循环步骤3.1.1.2，每移动一个时间间隔，计算一次相应的皮尔逊相关系数，直到为止。得到相关系数集合
[0123]
步骤3.1.1.4、从中找到最大值作为第k个缺陷区域像素点kxi与背景区域参考点温度热响应序列相关系数kρ
i&bg
：
[0124][0125]
步骤3.1.2、按照步骤3.1.1相同的计算方式，计算第k个缺陷区域像素点与缺陷特征实际区域参考点温度热响应的相关系数kρ
i&def
。根据相关系数kρ
i&bg
和kρ
i&def
得到第k个缺陷区域像素点xi的热扩散区域判定阈值kε
ttr_i
：
[0126][0127]
步骤3.2、根据计算得到的热扩散区域判定阈值kε
ttr_i
，当缺陷像素点与缺陷特征实际区域参考点之间的距离
kdwei
满足：
[0128][0129]
认为像素点xi位于热扩散区域，则去掉像素点xi。重复步骤3.2直到缺陷特征识别图像igc中第k个缺陷区域像素点xi都被访问到，去除所有位于热扩散区域的像素点xi，转至步骤3.3。
[0130]
步骤3.3、重复步骤3.1和步骤3.2直到k个缺陷区域像素点均被访问到，得到缺陷特征识别图像igc。
[0131]
步骤四、通过8-邻接freeman链码来检测缺陷特征识别图像igc中不同类型和不同位置缺陷的边界特征，对缺陷区域进行标定，实现缺陷像素粗定位。对第k，k＝1,
…
,k个缺陷区域进行像素标定定位的具体过程为：
[0132]
选择第k个缺陷区域的边界跟踪的起始像素点，记为ke0，以链码值dir＝0的方向为开始边界搜索的方向，逆时针方向旋转45
°
，按照方向dir＝1搜索下一个像素点。由于经过提取得到的缺陷特征图像为二值图像，在这个二值图像中，表征缺陷轮廓像素点的像素值
用“1”表示。基于此，在缺陷边界像素点的搜索过程中，如果像素点kx
t
的像素值为“1”，则认为像素点kx
t
是缺陷边界点，将其记为
ket
并保存。并将kx
t
的链码值赋值给它的前一个点像素kx
t-1
。将像素点kx
t
作为新的边界搜索起始点，逆时针方向旋转90
°
，按照方向dir＝2搜索下一个缺陷边界点。重复上述步骤，当搜索到起始点kb0时，缺陷区域搜索完毕。整个搜索过程中，检测到的边界点ke0,
…
,
kem
即为缺陷的边界轮廓。边界点集合kω
l
，kω
l
＝(ke0,
…
,
kem
)对缺陷进行了像素粗定位。重复上述过程直到k个缺陷区域访问完毕，得到k个缺陷区域的像素标定定位结果。
[0133]
步骤五、基于缺陷像素粗定位结果，根据提取到的缺陷边界像素点，对每一个标定的缺陷区域利用b
é
zier曲线拟合的方式实现亚像素缺陷特征定位。
[0134]
通过像素粗定位后，标定得到了k个缺陷区域，提取标定的第k，k＝1,
…
,k个缺陷区域的m个缺陷边界像素点，m个缺陷边界像素点构成像素点集合
kei
，i＝0,1,
…
,m。分别对相邻的缺陷边界像素点
kei
和像素点
kei+1
两两进行三次b
é
zier曲线拟合。由于i的最大值为m，因此像素点
kem+1
是不存在的。考虑到缺陷轮廓是封闭图形，用ke0作为
kem+1
与
kem
进行三次b
é
zier曲线拟合。利用三次b
é
zier曲线对缺陷边界像素点
kei
和
kei+1
两两进行拟合得到亚像素缺陷封闭轮廓的具体过程为：
[0135]
步骤5.1、在待拟合的边缘轮廓像素点
kei
和
kei+1
之间，构造控制点kki和
kci
。
[0136][0137][0138]
其中，
ka*
和
kb*
是任意给定的正数；(kxi,kyi)是缺陷边界像素点
kei
的位置参数。由于经过边界追踪算法检测得到缺陷的m个边界像素点
kei
(kxi,kyi)，i＝0,1,
…
,m，因此缺陷边界像素点集合中不存在像素点ke-1
、
kem+1
、
kem+2
，其对应的位置参数(kx-1
,ky-1
)，(kx
m+1
,ky
m+1
)，(kx
m+2
,ky
m+2
)也就不存在。考虑到缺陷轮廓是闭合的，用(kxm,kym)的值作为(kx-1
,ky-1
)的值，(kx0,ky0)的值作为(kx
m+1
,ky
m+1
)的值，(kx1,ky1)的值作为(kx
m+2
,ky
m+2
)的值；
[0139]
步骤5.2、根据b
é
zier曲线的定义，缺陷边界像素点
kei
和
kei+1
之间的三次b
é
zier曲线的表达式为：
[0140][0141]
根据三次b
é
zier曲线的表达式得到缺陷边界像素点
kei
和像素点
kei+1
之间的拟合曲线；
[0142]
步骤5.3、步骤四遍历缺陷边界像素点集合
kei
，i＝0,1,
…
,m所有的像素点，分别进行两两拟合得到相邻像素点之间的拟合曲线。当i＝m时，用ke0作为
kem+1
与
kem
进行三次b
é
zier曲线拟合；
[0143]
步骤5.4、所有缺陷边界像素点遍历结束，得到m条三次b
é
zier曲线。m条b
é
zier曲线平滑连接得到第k个缺陷区域的亚像素缺陷封闭轮廓kl。基于封闭轮廓kl实现缺陷亚像素
定位。
[0144]
步骤5.5、重复上述过程直到k个缺陷区域均被访问到，得到不同类型和不同位置的缺陷亚像素定位结果。
[0145]
步骤六、根据亚像素缺陷定位结果，对每一个标定的缺陷区域，对轮廓特征尺寸在水平方向和垂直方向，按照试件和图像的尺寸关系，分别进行映射以得到刻画缺陷不同特征参数的实际尺寸。对标定的第k，k＝1,
…
,k个缺陷区域，亚像素边界点的坐标记为(kxs,kys)，按照如下公式进行映射：
[0146][0147]
式中，l
×
w为检测试件的实际尺寸；m
×
n为拍摄得到的红外图像的尺寸，得到映射后的亚像素边界点坐标(kxs',kys')。
[0148]
步骤七、缺陷特征量化计算。周长和面积是常规且直接地反映缺陷的尺寸信息，基于几何结构特征参数数值，可直接量化评估各个缺陷特征区域的尺寸大小信息，进而体现各个缺陷特征区域之间的几何尺寸差异。缺陷的面积可以用来量化缺陷在图像中的大小。分别计算不同的缺陷区域的特征参数。对标定的第k，k＝1,
…
,k个缺陷区域，对将映射后的亚像素边界点(kxs',kys')集合勾画的缺陷实际尺寸中的边界轮廓记为kl。由于大多数缺陷特征区域呈现不规则性，因此缺陷边界轮廓kl的形状呈不规则图形。直接计算缺陷边界轮廓kl包围的缺陷区域的面积比较困难。对此，用缺陷边界轮廓kl包围的区域kω的像素点kpi的数目和作为面积特征参数。由于在步骤七中进行了边界尺寸映射，因此不再需要进行试件与图像之间的面积换算，此时统计得到缺陷边界轮廓kl包围的区域kω的像素点数目即为缺陷的面积特征参数：
[0149][0150]
通过统计缺陷轮廓包围的区域面积能够进一步确定缺陷的位置信息，达到标定缺陷的目的。在经映射后的缺陷轮廓图像中，基于缺陷边界轮廓特征进一步计算缺陷的周长、长径和短径特征。在计算缺陷的周长时，由于拟合得到的缺陷亚像素轮廓是通过按顺序两两拟合缺陷边界像素点
kei
和
kei+1
得到的，并且缺陷边界特征尺寸映射不会改变亚像素边界点之间的位置关系，因此通过曲线积分可以求得边界像素点
kei
和
kei+1
之间的拟合曲线的实际长度kpi：
[0151][0152]
其中，为缺陷边界像素点
kei
和
kei+1
之间的三次b
é
zier曲线的表达式。m个边界像素点两两之间的实际拟合曲线长度的和即为缺陷的周长kp：
[0153][0154]
同时由于缺陷特征区域的不规则性，基于亚像素缺陷边界轮廓得到的最小外接矩形框的尺寸同样也可以反映缺陷的形态特征。根据亚像素缺陷边界轮廓的拟合轮廓kl，找到其最小外接矩形，用最小外接矩形来刻画缺陷的分布情况。最小外接矩形的长边即为缺
陷的长径kdl，短边即为缺陷的短径kds。通过面积来利用等间隔旋转法来求解亚像素缺陷边界轮廓的最小外接矩形：亚像素缺陷边界轮廓图像在90
°
范围内等间隔地旋转，每次旋转都需要记录亚像素缺陷边界轮廓的外接矩形参数：分别找到亚像素缺陷边界轮廓在水平方向的坐标最大值kx
max
和最小值kx
min
，竖直方向的坐标最大值ky
max
和最小值ky
min
。则标定的第k个缺陷区域的外接矩形的长边为：
[0155]kdl＝max(ky
max-kymin
,kx
max-k
x
min
)
[0156]
短边为：
[0157]kds＝min(ky
max-kymin
,kx
max-k
x
min
)。
[0158]
则外接矩形的面积为kdl
×kds，本发明通过面积来定义最小外接矩形，计算最小外接矩形的具体过程为：首先根据亚像素缺陷边界轮廓分布坐标的最大值最小值能够得到一个没有旋转角度的亚像素缺陷边界轮廓区域的外接矩形，并计算当前状态下的最小外接矩形的长边、短边和面积，此时旋转角度α＝0
°
。设置等间隔旋转角度θ，对亚像素缺陷边界轮廓进行旋转一个角度θ，α＝α+θ，按照相同的方式求取亚像素缺陷边界轮廓的外接矩形，并记录相关参数。重复旋转过程，直到旋转角度α＝90
°
。最后，根据面积大小确定最小外接矩形。根据上述量化过程，能够基于亚像素缺陷边界轮廓拟合结果计算得到缺陷的物理特征参数。
[0159]
仿真验证：
[0160]
为了说明本发明的可行性和正确性，下面分别通过实验验证来对本发明的技术方案进行分析说明。
[0161]
本实施例中选择尺寸大小为200mm
×
200mm
×
5mm的高温陶瓷实验试件，材料亚表面嵌入了两种不同类型的圆柱损伤。左下角的是孔径φ为3mm的分层缺陷，其余圆柱均为填充了与高温陶瓷材料导热率不同的其他材料的夹杂缺陷。右上角的的是孔径φ为5mm，其余两个内部夹杂缺陷的孔径φ均为5mm。采集数据的红外相机分辨率为512
×
640，整个实验过程采集了99帧红外数据，红外数据以大小为512
×
640
×
99的数据块的形式进行存储。通过对红外热图像序列进行处理，得到表征缺陷的红外重构图像如图2所示。
[0162]
利用模糊c聚类分割的方式提取红外重构图像中的缺陷特征，缺陷特征提取图像如图3所示。
[0163]
在进一步定量识别缺陷特征的过程中，考虑到缺陷类型复杂和噪声的存在，缺陷边缘区域像素点可能存在错误划分的情况，更正确的缺陷区域划分结果能够提高识别的精度，这里基于缺陷特征提取图像进行了预处理，重新计算了边缘像素点的隶属度，根据隶属度对缺陷区域的边缘像素点重新进行了划分。后处理过程中，设置参数r＝3，n＝5，对缺陷特征提取图像进行后处理后得到缺陷特征检测图像如图4所示。
[0164]
根据后处理前后像素点的划分情况来看，后处理能够对缺陷边缘像素点进行重新划分，重新划分的结果与缺陷特征实际表征情况是相符的，说明了后处理方案的有效性。由于热扩散效应客观存在，通过特征提取得到的缺陷特征图像实际上包含了缺陷区域和热扩散区域，而后处理过程并不能去除热扩散区域的影响。考虑到热扩散程度未知，因此对位于缺陷区域的像素点的瞬态热响应曲线进行了对比分析。对缺陷区域进行标记，标记结果如图5所示，再分别去除每个缺陷区域的热扩散区域。对每个缺陷区域，将缺陷区域中心像素点作为缺陷特征实际区域参考点，计算其余像素点与参考点之间的距离，找到距离最大值dwemax
。通过计算各个缺陷区域像素点的瞬态热响应曲线与背景区域参考曲线的最大相似度ρ
i&bg
和与缺陷区域参考曲线的最大相似度ρ
i&def
来确定其阈值ε
ttr_i
。各个区域的部分阈值ε
ttr_i
计算结果和距离最大值d
wemax
计算结果如表1所示。
[0165]
表1各缺陷区域热扩散区域判定部分参数计算结果
[0166] ε
ttr_1
ε
ttr_2
ε
ttr_3
ε
ttr_4
ε
ttr_5
ε
ttr_6dwemax
缺陷区域a0.49930.49920.49920.49960.49940.499513.8996缺陷区域b0.50020.50010.50030.50010.49990.50003.1353缺陷区域c0.49990.49970.50000.50020.50020.500410.2249缺陷区域d0.50040.50020.50180.50040.50020.50033.0964
[0167]
基于热扩散区域判定参数去除经后处理的缺陷特征提取图像中的热扩散区域，得到图6所示的缺陷特征识别图像。经热扩散区域去除前后，各缺陷区域像素点数目对比如表2所示。
[0168]
表2热扩散去除前后各缺陷区域像素点数目对比
[0169] 缺陷区域a缺陷区域b缺陷区域c缺陷区域d去除热扩散区域前56786497131去除热扩散区域后1615916160
[0170]
根据表2的对比结果可以看出，本发明设计的通过时序排齐计算最大相关系数的方法确实能够达到修正缺陷区域的目的。
[0171]
提取缺陷特征识别图像的边缘轮廓信息，再进行边界跟踪。边界像素点构成的集合作为缺陷像素粗定位结果。选择a缺陷区域位于(123,199)的像素点、b缺陷区域位于(427,218)的像素点、c缺陷区域位于(428,421)的像素点和d缺陷区域位于(121,411)的像素点分别为各缺陷区域边界追踪起点，得到缺陷区域边界提取图像如图7所示。通过边界追踪提取到a缺陷区域44个边界像素点，b缺陷区域23个边界像素点，c缺陷区域39个边界像素点，d缺陷区域24个像素点。基于缺陷边缘轮廓像素点构成的集合对缺陷边缘轮廓进行拟合，设置拟合参数a＝0.25，b＝0.25，求得为拟合各缺陷区域轮廓而构造的部分控制点集如表3所示。
[0172]
表3各缺陷区域轮廓拟合部分控制点集
[0173][0174][0175]
得到亚像素缺陷区域定位结果如图8所示。从图8所示的亚像素缺陷区域定位图像来看，拟合出来的缺陷区域形状受缺陷边界像素点粗定位影响，拟合得到的缺陷轮廓连续光滑。放大后像素块出现失真的情况，而拟合曲线仍然保持连续，具有稳定的缺陷特征。
[0176]
进一步定量计算缺陷物理特征。首先进行尺寸映射，对大小200mm
×
200mm的试件进行检测，实验获得大小为640
×
512的红外重构图像。经过尺寸映射后，计算得到各项参数的像素尺寸即为缺陷的真实尺寸，高温陶瓷实验试件缺陷定量计算结果如表4所示。
[0177]
表4高温陶瓷实验试件缺陷定量计算结果
[0178][0179]
根据定量结果可以看出，本发明提出的缺陷提取和缺陷识别方法在识别缺陷特征方面具有一定的准确度。所有缺陷特征的误差范围在0.72％～6.75％，其中面积误差最小。说明了本发明设计的自动确定阈值的热扩散区域剔除方法能够根据缺陷特征不同、加热情况不同、热扩散程度不同自动调节阈值，得到的缺陷区域能够比较接近真实的缺陷区域。算法自动化程度提高的同时也提高了定量识别的准确度。
[0180]
对本发明的应用、修改和变化对本领域的技术人员来说是显而易见的。
[0181]
尽管本发明的实施方案已公开如上，但其并不仅仅限于说明书和实施方式中所列运用，它完全可以被适用于各种适合本发明的领域，对于熟悉本领域的人员而言，可容易地实现另外的修改，因此在不背离权利要求及等同范围所限定的一般概念下，本发明并不限于特定的细节和这里示出与描述的图例。