CN103308537B - 一种递变能量x射线成像图像融合方法 - Google Patents

Abstract

本发明公开了一种递变能量X射线成像图像融合方法，该方法从信息论角度出发，建立图像质量评价标准，并提取检测对象在各能量下X射线成像图像的有效区域；同时依据动态范围无约束下射线能量、图像灰度以及检测对象之间的唯一性，构建能谱图像序列之间的灰度变换模型，在保证图像灰度阶不会出现混乱的基础上实现多谱序列的有效区域的融合。该方法在完成成像系统动态范围扩展的基础上，降低了融合过程中对先验知识和人工参与的依赖性，提高了其智能性，对于推广变能量X射线成像技术在工程中的应用有重要意义。

Description

一种递变能量X射线成像图像融合方法

技术领域

本发明涉及X射线数字成像技术，特别是一种递变能量的X射线成像图像序列融合方法。

背景技术

在X射线成像过程中，扫描参数管电压和管电流决定了射线源的X光子通量密度和射线能量，直接影响辐射能量和射线图像质量；射线曝光量可以通过一个转换因子折算为物体的辐射剂量和吸收剂量。如果固定X射线成像系统的其他扫描参数，增加管电压可以增加辐射曝光量，增强射线的穿透力，同时穿过被检测物体的等效厚度也随之变长；增加管电流，可以提高X光子的通量密度，降低图像噪声，减少图像中的散斑块，提高图像质量。但是，由于数字化成像模式由于受射线转换效率以及光电转换效率的限制，即成像器件动态范围受限，所以对于有效厚度差异大的复杂结构件（结构复杂、长宽比较大、材料多组元化等），现有固定能量的射线成像模式无法在单一能量下实现检测对象结构信息的完整再现，易出现过曝光和欠曝光共存现象，导致投影信息严重缺失，从而无法进行X射线成像(digitalradiography,DR)检测。

针对上述问题，一种X射线成像装置和方法（CN：ZL200810161410.9）提出了一种递变能量X射线成像方式，即通过调节射线能量的方式，分别获取到检测对象的不同有效厚度分别对应的质量符合要求的射线图像，并通过对获取到的一系列射线图像进行图像融合，扩展图像的动态范围，实现了对有效厚度差异较大的检测对象的射线检测。

在现有递变能量成像方式中所涉及的图像融合主要有以下几种方法：一是依据大量的先验知识，寻找图像序列的最佳灰度带，将图像处于最佳灰度带的部分作为子图，依靠人工调节权系数对子图进行拼接获取融合图像，这种方法无法利用图像序列之间重叠区域，需要大量先验知识及人工参与，调节融合权系数，不适合工程应用；二是基于小波分析及类似算法的图像融合方法，由于这类方法未考虑DR图像序列中图像灰度逐渐增加的特殊性，会导致融合图像中灰度阶出现混乱，即被检测对象部分等效厚度大的区域的像素灰度值高于等效厚度小的区域的像素灰度值，与等效厚度越大，其灰度值越低的规律相悖；三是加权融合方法，这类方法可以保证融合图像中的灰度阶不出现混乱，但是如果约束条件或优化目标的设置不合理，则会造成融合权值的不合理，导致融合图像中出现“涟漪”状的伪边缘。

发明内容

本发明针对上述DR图像融合方法存在的缺点和不足，提供一种递变能量X射线DR图像融合方法，在完成成像系统动态范围扩展的基础上，降低了融合过程中对先验知识和人工参与的依赖性，提高了其智能性。

为了实现上述目的，本发明的技术方案是：

一种递变能量X射线成像图像融合方法采用递变能量X射线成像装置实现，递变能量X射线成像装置包括变剂量成像模块和图像处理模块，变剂量成像模块根据扫描过程中X射线透射方向上待检测物体有效厚度的变化，实时调整X射线源管电压对所述待检测物体进行扫描，然后将扫描得到的透射图像发送给图像处理模块；图像处理模块对接收到的各透射图像进行图像融合以及动态范围提升，得到最终所需透射图像，图1所示，图像处理模块中的图像融合方法包括以下步骤：

(1)系统初始化：包括硬件初始化和软件初始化，硬件初始化选择图像采集的初始能量；软件初始化定义一个与硬件初始化构建图像对应的等大的模板C₀和融合图像F₀；

(2)DR图像采集：设定步进电压，逐幅采集DR图像；

(3)提取DR图像的有效区域：结合灰度范围宽度及灰度直方图均衡性提取DR图像的有效区域；

(4)有效区域融合：依据图像间灰度变换模型进行有效区域的融合，判断融合图像信息是否完整，若已完整，停止采集，否则，继续采集图像进行融合；

(5)信息完整性判断：

获取融合图像F_j后，得C_j。若C_j=(0,0,…,0)，则结束图像采集，否则，继续调节电压，采集图像。

本发明具有的有益效果为：充分利用了图像序列之间的重叠区域，降低了融合过程中对先验知识和人工参与的依赖性，提高了图像融合质量及融合方法的智能性，同时可进一步为智能X射线成像的能量自适应调节提供技术支撑。

附图说明

图1为本发明递变能量X射线成像图像融合方法流程图。

具体实施方式

以下结合附图介绍本发明详细技术方案：

本方明是针对递变能量X射线成像装置中图像处理部分存在的不足提出的一种图像融合方法。递变能量X射线成像装置包括变剂量成像模块和图像处理模块，变剂量成像模块根据扫描过程中X射线透射方向上待检测物体有效厚度的变化，实时调整X射线源管电压对所述待检测物体进行扫描，然后将扫描得到的透射图像发送给图像处理模块；图像处理模块对接收到的各透射图像进行动态范围提升以及图像融合，得到最终所需透射图像。

本发明公开的递变能量DR图像融合方法的具体流程如图1所示，首先对系统参数初始化，然后采用合理的步进电压，逐幅采集DR图像，根据有效区域评判标准提取DR图像的有效区域，依据图像间灰度变换模型进行有效区域的融合，判断融合图像信息是否完整，若已完整，停止采集，否则，继续采集图像进行融合。

(1)系统初始化

系统初始化包括硬件初始化和软件初始化。硬件初始化主要是选择图像采集的初始能量，例如：通过电压、电流等。设构件图像为X=(x₀₁,x₀₂,…,x_0M)，m=1,2,…,M，在软件初始化中，定义一个与X对应的等大的模板C₀=(c₀₁,c₀₂,…,c_0M)=(1,1,…,1)和融合图像F₀=(f₀₁,f₀₂,…,f_0M)=(0,0,…,0)。

(2)DR图像采集

设定步进电压，逐幅采集DR图像。为了保证图像信息的完整性，该实例中电压步进值设置5KV。

(3)提取DR图像的有效区域

在各能量下DR图像中，由于射线能量与检测对象不匹配，在成像器件动态范围受限的情况下，易出现欠曝光和过曝光共存现象，只有有效区域能提供有效信息。因此，多谱序列有效区域的提取及质量直接关系到后续的多谱融合质量。

检测对象的结构信息最终以灰度形式反映在图像中，因此，对每幅图像可以选择一段范围内的灰度区域作为有效区域。考虑图像有效区域的性质：一方面，检测对象不同结构需要通过其对应像素间的灰度差来表现，这表明有效区域应该有一个较宽的灰度范围，能够表征一定的灰度差；另一方面，为提高图像质量，经常需要对其增强，其常用方法为灰度直方图均衡化，一定程度上图像的灰度直方图越均衡，其质量越高。图像信息熵方面，设图像总共有n个灰度阶，灰度阶由低到高为Hi(i=1,2,…,I)，灰度阶Hi在图像中出现概率为p_i，其中图像信息熵E可表示为由Lagrange乘数法可以得到时E取得最大值，当时，即意味着图像的灰度直方图最为均衡。由上述分析，可结合灰度范围宽度及灰度直方图均衡性设计图像的有效区域提取方法。

设图像中灰度范围为H_i1～H_i2区域，其宽度为H_i2-H_i1。记图像灰度直方图为S_i=S(H_i)，灰度阶H_i的均衡性可以用其δ邻域[H_i-δ,H_i+δ]内对应的S_i的标准差衡量，记为σ(H_i)（δ为一个正整数）。σ(H_i)越小，说明灰度阶H_i处均衡性越好。灰度直方图在上的平均均衡性可以定义为考虑有效区域的灰度范围宽度及其均衡性，可认为其满足优化模型

${\min \left(\frac{\underset{i=i_1}{\overset{i_2}{\mathrm{\Σ}}}\mathrm{\σ}\left(H_i\right)}{H_{i_2}-H_{i_1}}\right)}^{\mathrm{\α}}{\left(\frac{1}{H_{i_2}-H_{i_1}}\right)}^{\mathrm{\β}}$

s.t.S_i=S(H_i),（1）

$H_{i_2}-H_{i_1}>0,$

其中，α，β为正实数，其相对大小决定了优化目标中灰度直方图均衡性与灰度范围宽度两个因素的相对重要性，由于是非负数，α，β都是正数，令上述优化模型可进一步改写为：

$\min \frac{\underset{i=i_1}{\overset{i_2}{\mathrm{\Σ}}}\mathrm{\σ}\left(H_i\right)}{{(H_{i_2}-H_{i_1})}^{\mathrm{\γ}}}$

s.t.S_i=S(H_i),（2）

$H_{i_2}-H_{i_1}>0,$

求解式（2），得图像有效区域的灰度范围。

若σ(H_i)为连续函数，则取极值时有由于此处σ(H_i)是离散的，所以可遍历满足的进行求解，也可直接遍历图像中所有灰度阶求解。

为确保有效区域的高质量，结合人眼视觉系统对提取有效区域做质量评价，当有效区域达到一定的质量指标时，允许其进入下一步的融合过程。对比敏感度函数（ContrastSensitivityFunction，CSF）能较全面、客观地评价视觉功能，Mannos和Sakrison等人通过大量实验建立了CSF的近似曲线，其描述如下：

A(f)=2.6(0.0192+0.114f)exp(-(0.001f)^1.1)（3）

式中，空间频率（周期/度），其中，f_x、f_y分别为水平和垂直方向的空间频率。CSF具有带通滤波特性，对中间频段较为敏感。因此，采用小波分析方法对有效区域进行小波分解，利用其频域系数构建图像质量评价模型。有效区域是按照灰度范围选出的，这可能导致有效区域的不规则性。因此，随机选取有效区域内一定数量的图像块，对其进行评价，计算其平均质量，作为有效区域质量。

随机选取有效区域内N个P×P的图像块并对其归一化，记得到的归一化图像块为Z_n(n=1,2,…,N)，对Z_n进行L(L=4～5)级小波分解，得到L+1个频带，其中一个低频带，L个高频带，低频带中只有一个低频分量，每个高频带包含三个高频分量，记第l(l=1,2,…,L+1)个频带的第k(高频带中k=1,2,3，低频带中k=1)个分量为I_l,k(s,t)，s=1,2,…,S，t=1,2,…,T。计算各分量的方差var、梯度g和熵e，其中， ${\mathrm{var}}_{l,k}=\frac{1}{\mathrm{ST}}\underset{s=1}{\overset{S}{\mathrm{\Σ}}}\underset{t=1}{\overset{T}{\mathrm{\Σ}}}{(I_{l,k}(s,t)-\frac{1}{\mathrm{ST}}\underset{s=0}{\overset{S}{\mathrm{\Σ}}}\underset{t=0}{\overset{T}{\mathrm{\Σ}}}{I}_{l,k}(s,t))}^2,$

$g_{l,k}=\frac{1}{(S-1)(T-1)}\underset{s=1}{\overset{S-1}{\mathrm{\Σ}}}\underset{t=1}{\overset{T-1}{\mathrm{\Σ}}}\sqrt{{(I_{l,k}(s,t)-I_{l,k}(s+1,t))}^2+{(I_{l,k}(s,t)-I_{l,k}(s,t+1))}^2},$

然后将各频带内所包含分量的计算结果分别求和，对高频分量， ${\mathrm{Var}}_l=1n\left(\frac{1}{3}\underset{k=1}{\overset{3}{\mathrm{\Σ}}}{\mathrm{Var}}_{l,k}\right),$ $g_l=1n\left(\frac{1}{3}\underset{k=1}{\overset{3}{\mathrm{\Σ}}}{g}_{l,k}\right),$ $e_l=1n\left(\frac{1}{3}\underset{k=1}{\overset{3}{\mathrm{\Σ}}}{e}_{l,k}\right),$ 对低频分量，Var_l=ln(var_l,k)，g_l=ln(g_l,k)，e_l=ln(e_l,k)得到Z_n在第l个频带的综合评价指标根据归一化的CSF选择合适的权值vl，得In的质量进一步得到有效区域质量选择合适的阈值Q'，当Q>Q'时，认为有效区域可参与融合过程，否则，认为有效区域无效，继续采集图像。

(3)有效区域融合

融合图像可等效为在没有受到成像系统灰度范围限制的情形下被检测对象在某一电压下的DR图像。对于同一被检测对象，由于其信息的同一性，当成像系统的灰度范围足够时，其在不同电压下的DR图像之间必定存在一个变换关系。可通过不同电压下DR图像有效区域的重叠区域的灰度值变化寻求这一变换关系。依据这一变换关系，将不同电压下DR图像有效区域的灰度值变换成同一电压下的灰度值，然后对有效区域拼接，即得到被检测对象在对应电压下处于相应区域的等效DR图像。

设电压U_j(j=1,2,…,J)下DR图像X_j=(x_j1,x_j2,…,x_jM)中第j₁个像素，第j₂个像素，…，第个像素，共M_j(j=1,2,…,J)个像素构成其有效区域，记这些像素的位置集合为B_j，依据位置集合B_j可以得到图像Y_j=(y_j1,y_j2,…,y_jM)，公式为 $y_{\mathrm{jm}}=\left\{\begin{array}{cc}{x}_{\mathrm{jm}},& m\∈B_j\\ 0,& m\∉B_j\end{array}\right.,$ (m=1,2,…,M)，有效区域融合后图像记为F_j=(f_j1,f_j2,…,f_jM)。当被检测对象中存在孔洞时，其对应区域在最低电压下的图像X₁中表现为过曝光区域，为了表示统一，将此过曝光区域的像素对应位置集合记为B₀，类似于有效区域的位置集合，依据B₀更新F₀，公式为 $f_{0m}=\left\{\begin{array}{cc}{x}_{1m},& m\∈B_0\\ 0,& m\∉B_0\end{array}\right.,$ (m=1,2,…,M)，继而计算F₁，F₁=F₀+Y₁。

设r∈B_j-1∩B_j={1,2,…,R}(j>1)，对y_j-1,r和y_j,r进行数据分析发现有关系式y_j,r≈α_j-1,jy_j-1,r+b_j-1,j成立，其中α_j-1,j，b_j-1,j为参数，且α_j-1,j>0。记(y_j-1,r,y_j,r)到直线y=a_j-1,jx+b_j-1,j距离为d_r,j-1,j，考虑到Y_j-1与Y_j都存在噪声，故选择a_j-1,j，b_j-1,j，使其满足 $\underset{r=1}{\overset{R}{\mathrm{\Σ}}}{\left(d_{r,j-1,j}\right)}^2$ 最小，即 $\min \underset{r=1}{\overset{R}{\mathrm{\Σ}}}\frac{{(a_{j-1,j}{y}_{j-1,r}-y_{j,r}+b_{j-1,j})}^2}{1+{\left(a_{j-1,j}\right)}^2}.$

解得a_j-1,j，b_j-1,j后，融合F_j-1与Y_j。记A₀=B₀，令A_j=A_j-1∪B_j(j>0)。F_j-1与Y_j融合后的图像F_j满足 $f_{\mathrm{jm}}=\left\{\begin{array}{cc}{a}_{j-1,j}{f}_{j-1,m}+b_{j-1,j},& m\∈A_{j-1}\\ y_{{\mathrm{jm}}^{,}}& m\∈(B_j-A_{j-1})\end{array}\right.,$ (m=1,2,L,M)。

(4)信息完整性判断

获取F_j后，计算c_jm(m=1,2,…,M)，计算公式为c_jm=c_0mg(f_jm)， $\left(x\right)=\left\{\begin{array}{cc}0,& x>0\\ 1,& x=0\end{array}\right.,$ 由c_jm构成C_j=(c_j1,c_j2,…,c_jM)。若C_j=(0,0,…,0)，图像采集结束，否则，继续调节电压，采集图像。

Claims (1)

1.一种递变能量X射线成像图像融合方法，采用递变能量X射线成像装置，所述递变能量X射线成像装置包括变剂量成像模块和图像处理模块，其中：变剂量成像模块根据扫描过程中X射线透射方向上待检测物体有效厚度的变化，实时调整X射线源管电压对所述待检测物体进行扫描，然后将扫描得到的透射图像发送给图像处理模块；图像处理模块对接收到的各透射图像进行图像融合以及动态范围提升，得到最终所需透射图像，其特征在于，图像处理模块中的图像融合方法包括以下步骤：

(1)系统初始化：包括硬件初始化和软件初始化，硬件初始化选择图像采集的初始能量；软件初始化定义一个与硬件初始化构建图像对应的等大的灰度为1的模板C₀，以及灰度为0的融合图像F₀，C₀＝(c_0l，c₀₂,...,c_0M)＝(1，1,…，1)，F₀＝(f_0l，f₀₂,...,f_0M)＝(0，0,…，0)；

(2)X射线成像图像采集：设定步进电压，逐幅采集X射线成像图像；

(3)提取X射线成像图像的有效区域：结合灰度范围宽度及灰度直方图均衡性提取X射线成像图像的有效区域；

(4)有效区域融合：依据图像间灰度变换模型进行有效区域的融合，判断融合图像信息是否完整，若已完整，停止采集，否则，继续采集图像进行融合；

(5)信息完整性判断：获取融合图像F_j后，计算c_jm，m＝l,2，...，Μ，计算公式为c_jm＝c_0mg(f_jm), $g\left(x\right)=\{\begin{array}{cc}0,& x>0\\ 1,& x=0\end{array},$ 其中，C₀＝(c_0l，c₀₂,...,c_0M)，c_0m为初始模板C₀的像元，有效区域融合后图像记为F_j＝(f_j1,f_j2,…,f_jM)，f_jm为有效区域融合图像F_j的像元，g(f_jm)为融合图像F_j的像元标记函数，由c_jm构成C_j＝(c_jl，c_j2,...,c_jM)，若C_j＝(0，0,…，0)，则结束图像采集，否则，继续调节电压，采集图像。