郭 牛, 桂志国, 刘 祎, 康嘉琦
(1.中北大学 省部共建动态测试技术国家重点实验室, 山西 太原 030051;2.中北大学 信息与通信工程学院, 山西 太原 030051)
目前, 工业领域中一般使用X射线平板探测成像系统, 通过配套软件将生成的16位图像显示在8位PC显示器上, 图像内的大量信息在灰度范围从16位线性压缩到8位时都丢失了[1]。 此外,由于X成像的对比度较低, 图像在黑暗区域、 高亮区域以及物体重合区域显示不清, 增加了缺陷检测的难度。 最常见的手动调窗技术费时费力且效果不好, 而使用图像处理技术可以实现自动化且快速提升图像整体的观感和特征区域的辨识度,使工件的轮廓与细节更加清晰, 达到真正的检测和探伤目的[2]。 因此, 增强图像有用信息的同时抑制其无用信息是工业射线图像处理中待解决的关键问题之一。
图像的对比度拉伸算法可以分为全局(空间不变)算法和局部(空间变化)算法, 它们的区别在于部分算法只基于图像全局特性进行调整, 另一部分算法结合图像的全局和局部的视觉特性进行不同的调整。 全局对比度拉伸算法对整幅图像的所有图像像素灰度值采用相同的变换函数, 也就是说整幅图像中的所有相同的像素灰度值, 在变换后的灰度值也是相同的。 但是由于统一不变的变换曲线不能自适应图像的不同区域, 会导致图像处理结果在细节、 亮度以及颜色上的损失。 局部对比度拉伸算法在图像中的不同空间区域采用不同的变换函数, 换言之, 拉伸前后图像像素灰度值相同的区域, 在拉伸后的灰度值可能不同, 这取决于其所在位置和周围的像素灰度值[3]。 局部算法更多得考虑人眼的视觉特性, 保留了图像中各个区域特有的信息并带来了较好的视觉效果,但是其算法计算量较大, 同时还容易给图像带来噪声、 光晕及灰度值跳变, 提高了算法的复杂性。为了解决这些问题, 学者们提出了各种类型的局部对比度增强算法。2011年,Duan等[4]提出了基于局部对比度拉伸的色调映射算法, 该方法使图像具有更加丰富的层次和细节, 通过加入双边滤波处理, 很好地消除了光晕现象和噪声影响;2018年, 徐寒等[5]提出了基于小波变换的工业X射线图像模糊增强算法, 提出了新的隶属度函数,较传统模糊增强算法加大了变换强度, 提升了运算速度;2019年, 周冲等[6]提出了一种基于梯度场的工业图像增强方法, 该方法有效地增强了图像细节信息, 使平均梯度和清晰度都有所提升;2020年, 许纪亚等[7]提出了将拉普拉斯金字塔增强与二维直方图均衡化(2DHE)算法相结合, 在保留2DHE增强对比度效果的同时, 弥补了细节丢失的不足;2021年,Liu等[8]将梯度域色调映射和模糊逻辑相结合应用于工业X射线图像增强, 此方法对复杂工件X射线图像的缺陷检测具有良好的视觉效果。
受上述研究启发, 本文提出了一种基于模糊熵的X 射线图像局部对比度拉伸增强算法, 在16位X射线图像处理上, 本文算法能够展现更多的图像信息和细节, 并在一定程度上抑制均匀区域的噪声, 同时能够很好地解决边缘区域的某些过冲现象, 提高图像的整体质量, 实验结果证明了本文所改进算法的有效性。
Duan等提出了基于局部对比度拉伸的色调映射算法, 该算法的基本原理是: 首先对高动态范围图像进行对数映射来压缩动态范围, 然后进行局部对比度拉伸, 算法中为了消除光晕现象和噪声影响, 使用了双边滤波处理, 该算法结果呈现出较好的视觉效果, 算法核心部分是所提出的局部对比度的拉伸[4]。
该算法可由下式表示
式(1)是对高动态范围图像进行对数映射操作, 其中:Lmin和Lmax分别为原始图像素灰度值的最小值和最大值,P(I)为原始图像,t是为了保护对数函数, 其取值范围为[10-6,1],D(I)为对数映射图像。 在式(2)中:D(I)min为局部区域像素灰度范围拉伸前的最小值,D(I)max为局部区域像素灰度范围拉伸前的最大值,D(I)in为局部区域像素灰度范围拉伸前的中心像素灰度值,D(I)'
min和D(I)'max分别为局部区域像素灰度范围拉伸后的最小值和最大值,D(I)out为局部区域像素灰度范围拉伸后的中心像素灰度值。 在式(3)和(4)中:D(I)ming和D(I)maxg分别为原图经过对数映射后整幅图像的灰度最小值和最大值,β用于控制拉伸比例系数。在式(5)中,c为常数,取1.5。
该算法的工作原理可用图1来描述。 因人类的视觉系统对亮度的响应具有对数非线性, 所以, 该算法先对高动态图像的亮度作全局对数变换, 将亮度范围归一化。 之后, 在图像每个像素点所在的局部块, 将其对比度(局部区域内亮度值最大的像素点和亮度值最小的像素点的差值)进行拉伸。 在图1中, 为了扩大每个窗口的对比度, 由式(3)和式(4),将D(I)min的值向D(I)ming方向移动, 将D(I)max的值向D(I)maxg方向移动。 式(5)中计算得到的β控制每个窗口对比度扩大的比例。
图1 基于局部对比度拉伸的色调映射算法工作原理图Fig.1 Schematic of the working principle of local contrast stretch based tone mapping
针对上述算法在16位X 射线图像上表现出的图像对比度较低、 图像细节不清晰、 边缘区域过冲等现象, 本文首先通过全局对数映射和双边滤波对图像进行预处理, 再引入模糊熵来自适应调节局部对比度的拉伸程度, 对图像不同区域的细节信息作出改善, 达到图像增强的目的。
本文算法流程图如图2所示。
图2 本文算法流程图Fig.2 Flowchart of the algorithm in this paper
1.2.1 预处理
由于在低Kvp(X射线能量)条件下获取的X射线图像中的有效信息通常是全暗的, 因此, 在图像对比度拉伸之前需要进行初步的全局映射。应用对数公式把高动态范围光谱图的亮度值映射在0~1之间, 对数映射公式如下
双边滤波是一种非线性滤波器, 它有保持边缘、 降噪平滑的效果。 同其他滤波原理一样, 双边滤波也是采用加权平均的方法, 用周边像素亮度值的加权平均代表某个像素的强度, 所用的加权平均基于高斯分布。 双边滤波的权重在考虑像素欧氏距离的同时考虑了邻域内像素灰度值的辐射差异, 在计算中心像素时同时考虑这两个权重,双边滤波的公式如下。
1.2.2 基于模糊熵的局部对比度拉伸
本文通过引入模糊熵自适应调节局部对比度的拉伸程度来展现图像细节信息和达到图像增强的目的。
文中选取使用3×3的窗口, 找到窗口中的最小值F(I)min和最大值F(I)max, 令窗口内的最大值与最小值的差为窗口的对比度, 通过调整窗口中的最小值和最大值达到拉伸对比度的效果, 工作原理图如图3所示。
图3 基于模糊熵的X射线图像局部对比度拉伸增强算法工作原理图Fig.3 Working principle diagram of local contrast stretch enhancement algorithm for X-ray images based on fuzzy entropy
式中:Q(I)min与Q(I)max分别为局部最小与最大像素灰度值拉伸后的最小值和最大值; 参数b是为了防止图像的对比度拉伸过度, 在待处理的图像中产生噪声; 参数γ控制图像的细节对比度拉伸程度, 使图像的细节位置对比度的拉伸程度随该位置的模糊熵变化来更贴近人类的视觉。2个参数的数学表达式如下
WANG 等[9]提出的图像模糊熵反映了窗口中心(i,j)邻域中像素灰度空间变化率的大小, 即当3×3窗口处于较为均匀区域时, 窗口中心点的模糊熵较小; 而当3×3窗口处于灰度值变化较大区域时, 窗口中心点的模糊熵较大。 图像模糊熵公式具体如下式中:x(i,j)为对数映射图像坐标(i,j)处像素的灰度值;um(x(i+k,j+l))为隶属度函数, 表示在3×3窗口中像素值相对于窗口中心(i,j)处像素灰度值的隶属程度。 窗口Wn(i,j)中(i+k,j+l)处的像素灰度值与(i,j)处的像素灰度值的差越小,(i+k,j+l)处像素灰度值的隶属度越大; 反之, 该像素灰度值的隶属度越小。 隶属度函数公式具体如下
由式(14)和图4得知,16位X射线图像中,γ值小于1, 进一步得知当图像像素点处于均匀区域, 即图像模糊熵eA较小时,γ较大; 当图像像素点处于灰度值变化较大区域, 即图像模糊熵eA较大时,γ较小。 由式(11)与式(12)得知, 当γ较大时, 即图像均匀区域对比度拉伸较小; 当γ较小时, 对图像灰度值低且变化较大区域的对比度拉伸较大, 从而在展现出更多图像细节的同时,抑制了均匀区域的噪声。 式(14)中不同L值下图像处理结果(以16位X射线图像1为例, 因为该图像细节信息更为复杂, 增强难度较大)如图5(a)~图5(f) 所示。
图4 不同L 值下的γ 曲线Fig.4 γcurve at different L values
图5 16位X射线图像1的处理结果Fig.5 Processing result of 16-bit X-ray image 1
式中:H(I)为全局色调映射图的像素灰度值,F(I)min与F(I)max分别为滤波图中局部像素灰度值的最小值与最大值,Q(I)min与Q(I)max分别为局部最小与最大像素灰度值拉伸后的最小值和最大值,Q(I)为输出结果图的像素灰度值。
本文采用的编程环境为Visual Studio 2019+opencv4.0.0, 在硬件配置Intel(R)Core(TM)i7-10700 CPU@2.90 GHz,16 G 内存的计算机上完成实验仿真。
NRSS、SMD2和EOG 的值越高, 图像的细节信息越丰富, 信息熵Entropy越高表明图像包含信息越多。
通过设置式(14)中的参数L观察其对图像对比度拉伸的影响, 分别选取L=0.2,0.3,0.35,0.4,0.5, 其实验结果图分别如图5(a)~图5(f)所示, 图5(a)~图5(f)的客观评价指标如表1所示。
表1 16位X射线图像1在本文方法中取不同L 值的图像客观评价指标Tab.1 Objective evaluation index of 16-bit X-ray image 1 when taking different L values in the proposed method
对图5(a)~图5(f)的主观视觉分析明显看出,通过调整参数L的大小,16位X射线图像1随参数L的增大, 图像亮度更强。 再对图5(a)~图5(f)分别进行客观评价指标分析, 在表1中用粗体突出显示16位X射线图像1使用本文方法时, 选取不同L值的图像质量指标最佳值, 可以明显看出当L=0.35时, 图像的NRSS,SMD2,EOG,Entropy 4项质量指标均大于其他L值效果, 为图像质量指标最佳值。 结合主观视觉分析与客观评价指标分析, 当参数L超过0.35时, 随着L的增大, 图像亮度提升的同时, 图像客观评价指标降低, 因此, 本文将L=0.35作为参数L的上限。
为了验证本文算法的有效性和实用性, 将其与自适应直方图均衡(AHE)和文献[4]方法结果进行比较, 将算法应用到多种图像上进行实验分析, 在主观视觉效果和客观质量评估方面评价算法的优劣。
2.2.1 主观视觉效果分析
图6~图8分别展示了对3种16位X射线图像的处理结果。
图6 16位X射线图像1的处理结果Fig.6 Processing result of 16-bit X-ray image 1
图8 16位X射线图像3的处理结果Fig.8 Processing result of 16-bit X-ray image 3
为了更好地观察效果, 将感兴趣区域(ROI)(例如缺陷处)放大显示在第2行, 在原始X射线图像中,ROI区域用矩形框标记。 明显看出原始图像ROI区域内, 几乎看不见任何细节, 使用AHE方法和文献[4]方法得到的结果图中, 图像整体亮度得到明显提升, 但是缺点同样很明显:(1)图像处理后的整体视觉效果,AHE方法出现了区块效应(如图6(b)、 图7(b)和图8(b)第1行的箭头处); (2)AHE方法与文献[4]方法处理结果的对比度较低, 不能够很好地观察到原始图像中存在的缺陷现象(如图6(b)、 图6(c)、 图7(c)、图8(b)和图8(c)第2行的椭圆区域内); (3)AHE方法的边缘区域过渡不明显(如图7(b)第2行的椭圆区域内), 且在图像均匀区域产生了较多的噪声(如图7(b)第2行的箭头处); (4) 文献[4]方法图像边缘过渡区域存在一定的过冲现象(如图8(c)第2行的箭头处)。
图7 16位X射线图像2的处理结果Fig.7 Processing result of 16-bit X-ray image 2
就主观视觉效果而言, 本文方法相较于AHE方法与文献[4]方法, 在结果图上具有较好的对比度拉伸效果, 在图像均匀区域抑制了噪声现象,且一定程度上削弱了边缘过渡区域的过冲现象。
2.2.2 客观评价指标
在表2中, 同样使用粗体突出显示实验方法的质量指标最佳值, 可以明显看出, 本文方法对16位X 射线图像1(全)、16 位X 射线图 像1(ROI)、16位X射线图像3(全)、16位X射线图像3(ROI)处理后, 增强结果图的客观评价指标均为最佳值, 而对于16位X射线图像2(全)、16位X射线图像2(ROI), 通过结合图7(b)的第1行和第2行明显看出,AHE方法在图像中产生较多噪声与区块效应, 使得客观评价参数计算过大, 故对于16 位X 射线图像2(全)和16 位X 射线图像2(ROI)在不同方法下的增强中,AHE方法的客观评价指标为最佳值。
表2 X射线图像客观评价指标Tab.2 Objective evaluation index of X-ray images
通过主观视觉效果分析与客观评价指标数据的对比, 得知本文算法较AHE方法和文献[4]方法在16位X 射线图像上的增强效果取得了较好的清晰度和更客观的评价。
本文提出了一种改进的局部对比度拉伸算法,该算法通过图像模糊熵来控制局部对比度的拉伸。从经过处理的16位X射线图像效果可以看出, 与文献[4]提出的局部对比度拉伸算法相比, 本文算法的增强效果图中噪声更少, 且一定程度上减弱了边缘过冲现象, 增强了图像的真实感; 图像信息更多, 细节更加明显, 图像的整体质量得到了明显地提升。 本文方法目前对滤波方法存在改进空间, 对于参数的选取未能实现自适应确定。 在未来考虑引入L0梯度最小平滑滤波进行尝试, 结合图像局部信息尝试实现参数自适应操作。