高速飞行条件下长曝光图像复原方法

冯雨欣,厉小润,丁楫刚

(浙江大学电气工程学院,浙江 杭州 310027)

1 引 言

随着精确制导需求的发展,光学成像探测跟踪技术已成为必然趋势[1]。而在大气层、临界层中高速飞行的飞行器,其周边光学成像环境复杂,飞行器目标图像在获取和传输的过程中受到多种因素影响,包括变化流场引起的气动光学效应、成像设备光学系统导致的离焦模糊和衍射、成像设备与目标之间的相对位移、图像传输噪声、热辐射噪声等,这些因素会导致图像质量的退化,进而引起图像的模糊失真与附加噪声,严重降低图像的信噪比和目标的显著度,影响后续的目标探测精度[2]。因此,对退化图像进行去噪、去模糊处理具有重要的研究意义。

图像复原首先需要考虑噪声抑制问题。部分研究在建立噪声先验模型后,添加正则化项抑制噪声。文献[3]构造低秩去噪模型,利用奇异值分解求解模型;文献[4]构造非凸变分模型,交替最小化求解清晰图像;文献[5]结合全变差正则化罚项抑制加性噪声。但该类方法通常假设噪声符合高斯加性噪声模型,无法有效去除非高斯模型噪声及残留噪声;同时,去噪结果依赖正则化项的选择,目前仍未有通用的正则化参数选择准则,部分算法仍依赖手工调参。基于此,部分研究将去噪作为单独的预处理步骤,最大程度减少噪声对后续去模糊复原的影响。文献[6]提出改进加权核范数最小化方法,结合噪声强度自适应调整参数;文献[7]提出三维块匹配滤波算法,结合自适应基底,充分利用图像非局部特征,同时引入导向滤波。形态学去噪由于其计算复杂度低、可扩展性高的特点,在单一背景的图像中取得了良好的去噪效果[8-10],但单一固定的形态学结构元素会损失图像信息,造成噪声残留。

可使用点扩散函数(PSF)描述图像的模糊退化过程[11]。去模糊方法可分为两类,一类方法是求解PSF的同时迭代估计清晰图像的方法,该方法适用于PSF先验知识较少的情况,如文献[12] 通过编码后的模糊算子去逼近一个未知的模糊算子并搜索对应的清晰图像,文献[13]提出的基于信息熵的迭代盲反卷积算法;另一类方法是首先求解PSF,再进行反卷积估计清晰图像,该类方法适用于PSF有先验模型或较多先验知识的情况,如文献[14]提出基于稀疏先验模型的复原框架,文献[15]提出的基于块强度分布算法,文献[16]提出的自适应盲复原算法,该方法使用模糊图像频谱和能量约束求解PSF,但由于实际图像频谱仅在一定范围逼近对称,尤其容易受到运动模糊条纹干扰[17],使用固定阈值或常数代替部分频谱信息将无法充分利用已知先验,造成PSF估计误差。

由于高速飞行条件下长曝光图像可以近似认为部分PSF符合G类点扩散函数模型,因此选择先求解PSF后求解清晰图像的方法。针对上述现有复原算法存在的各种问题,本文提出一种面向长曝光图像的高速飞行条件下图像复原方法,首先研究高速飞行条件下图像退化机理,提出图像复原框架;构造多尺度多方向形态学滤波器,计算噪声特征图像计算结构元素尺度;融合二维小波的多分辨率特性,对低频子频段完成去模糊后再与高频子频段小波重建,进一步抑制残留噪声;使用结合最大类间方差法的Radon变换计算运动模糊点扩散函数参数,去除运动模糊伪影;最后构造多直线拟合模型,选取图像频谱的多方向直线拟合,同时计算基于图像梯度的对称频谱估计,进而求解G类点扩散函数参数,完成图像复原。

2 高速飞行条件下图像退化机理

图像复原是图像退化的逆过程,因此研究图像退化机理尤为重要。可以使用点扩散函数h(x,y)与噪声n(x,y)描述图像退化过程[18]:

g(x,y)=f(x,y)⊗h(x,y)+n(x,y)

(1)

G(u,v)=F(u,v)H(u,v)+N(u,v)

(2)

式中,f(x,y)、H(u,v)为清晰图像及其频谱,g(x,y)、G(u,v)为退化图像及其频谱。

在大气层、临界层中,高速飞行条件下的目标会首先受到气动光学效应的影响,产生图像模糊与畸变;目标与流场摩擦产生的气动热辐射噪声大幅降低图像信噪比,造成目标与背景的对比度下降;长曝光图像会记录探测器成像系统与目标之间的相对位移,造成图像产生某方向的伪影;成像系统自身成像机理导致的模糊和噪声也会造成图像退化。其图像退化机理框架如图1所示。

图1 高速飞行条件下图像退化机理框架图

3 高速飞行条件下图像复原

图2为图像复原流程图。由于图像受到多种因素的退化影响,且去模糊的关键在于PSF的计算,使用单个PSF表示多种模糊退化会增加PSF建模难度,本文对三部分退化分开计算,充分利用图像频谱的先验知识。

图2 高速飞行条件下图像复原流程图

3.1 自适应加权形态学去噪

形态学计算是定义在两个集合之间的集合论计算,其基本运算包含腐蚀、膨胀、开运算和闭运算。单一开-闭滤波器或闭-开滤波器会造成噪声脉冲的残留,使得形态滤波器的输出统计偏倚更严重;选择单一的结构元素会模糊图像细节,造成其余方向图像信息的丢失。本文采取多尺度多方向结构元素,构造自适应加权形态滤波器,其步骤可总结如下:

1)结构元素形状的选取。常用的方形和圆形结构元素会对图像边缘造成较大的模糊,因此本文选择不同方向不同尺寸的线性结构元素:

sp,q={s1,1,s1,2,…,s1,q,…,sp,1,…,sp,q}

(3)

其中,p为方向;q为尺寸。

2)结构元素方向q的选取。选取4个方向,分别为θ=0°、45°、90°、135°。

3)结构元素尺寸q的选取。尺寸为q的结构元素会滤除比q小的噪声而保留比q大的图像信息,因此准确计算噪声尺寸尤为重要。本文使用噪声特征图像确定噪声尺度。定义多尺度元素:

(4)

则尺度q对应的噪声特征图像为:

(5)

(6)

计算每个尺度特征图像对应的像素值,得到尺度分布:

(7)

其中,Ι(q)为尺度q对应的像素数目。根据尺度分布,确定噪声与图像信息的最佳尺度分界阈值TN,表示尺度q

q={1,2…,T}

(8)

对不同方向结构算子滤波结果加权处理,权值ωk计算方式[9]如下:

(9)

(10)

其中,f(x,y)为退化图像。

3.2 基于预实验的成像系统模糊PSF重建

探测器成像系统在成像过程中依赖光学部件之间的相互作用,其成像性能与质量对后续目标图像有直接影响,可以通过预实验进行评估。在成像系统可测条件下,以成像脉冲(小亮点)δ作为输入,则成像系统输出为:

Gc=HcA

(11)

其中,A为脉冲δ的强度;Gc为成像系统输出,可以求得成像系统PSF的频域表示:

(12)

由于实际实验中成像系统带有噪声,可固定机位多次拍摄,取平均后再进行处理。

3.3 基于改进Radon变换的运动模糊PSF重建

考虑长曝光图像运动符合单一方向线性匀速运动,其光学传递函数幅值包含sinc函数,模糊图像频谱图会呈现明暗相间的平行条纹[17],且:

(13)

(14)

其中,L为运动模糊PSF的尺度;θ为角度;φ为条纹倾角;d为条纹间距。因此,根据平行条纹可计算出运动模糊PSF的参数。

实际的图像频谱往往存在大量毛刺、噪点和中心十字亮线干扰,为解决上述问题,本文对采用基于图像质心的分割加权方法,以图像质心为中心,将图像分为4幅子图像,分别求取子图像的频谱图;将频谱图动态范围压缩、居中;使用最大类间方差法设定阈值,将灰度图像转换为二值图像,提取条纹结构,如图3所示。

图3 二值化图像频谱示意图

对图3(b)进行1°～180°的Radon变换,得到不同方向直线积分的投影值矩阵R,其最大值Rmax对应的角度即为条纹倾角φ,由式求解角度θ。将图3(a)顺时针旋转θ角度,计算各列像素的和Gc_col=∑colGc,其形状近似sinc函数,采用启发式算法[19]计算主瓣两侧极小值之间的距离2d,由式求解尺度L。求取不同子图的运动模糊PSF后求平均,则运动模糊点扩散函数可描述为:

(15)

3.4 气动光学效应PSF重建

长曝光退化图像的气动光学效应PSF符合G类PSF模型[20],在频域可表示为:

|G(u,v)|=exp(-α(u2+v2)β)

(16)

将式(16)代入式(2)中,并取模归一化得到:

-α(u2+v2)β=ln|G(u,v)|-ln|F(u,v)|

(17)

由于G类扩展函数呈辐射对称性,平面上过原点的任一直线都应满足。但在实际的退化图像中,由于噪声及其他随机因素的影响,频谱信息仅在一定程度上逼近辐射对称,不同直线拟合得到不同的结果,导致PSF估计失真。为解决上述问题,本文使用多方向直线进行拟合,得到多组PSF数值后平均,同时考虑到运动模糊图像频谱的频谱特征,选择u=0、v=0、v/u=tan(φ-π/2)三条直线,构造多直线平均拟合模型,提升参数求解的准确性。

同时,公式中清晰图像频谱F(u,v)为未知参数,文献[16]使用常数A代替,适应性差,无法充分利用退化图像频谱信息。本文通过对大量清晰图像频谱研究,提出基于图像梯度的对称频谱估计重建F(u,v),相较于先验固定参数的方法,本文方法具有自适应性。如图4所示,模糊图像与清晰图像的频谱图像在对称中心一定范围内较为相近,距离对称中心越远,频谱下降程度越大。因此,本文提出三段式拟合清晰图像频谱,其步骤可总结如下:

图4 清晰图像与模糊图像频谱对比

1)以对称轴中心作为原点建立坐标系。选取u=0直线方向上退化图像的频谱信息,得到G(v),且v∈(-N,N);

2)自适应确定阈值,T1取|G(v)|中心点附近第一个极小值点对应的v,该部分表示未受到点扩散函数影响的频谱;T2取图像梯度Δ|G(v)|极小值点对应的v;a取T2到N内|G(v)|下降的斜率,即:

(18)

则清晰图像频谱直线可表示为:

(19)

3)根据式计算-α(v2)β,仅保留原点两侧首个极小值点间的点,采用二次拟合计算参数α和β。

4)重复上述步骤,得到v=0、v/u=tan(φ-π/2)方向上的参数,求和平均后计算,得到气动光学效应的PSF。

3.5 自适应参数滤波

1)一次维纳滤波,使用参数ξ代替功率谱参数估计:

(20)

(21)

由于图像在局部小尺度上较为平稳,因此采用滑动窗口计算局部方差:

(22)

其中,PQ为局部小尺度窗口区域;μ(x,y)为窗口局部均值。

2)计算噪声方差:

(23)

其中,f′(x,y)表示一次维纳滤波的复原结果。

3)计算f′(x,y)的功率谱,带入:

(24)

4 仿真结果与分析

为验证本文算法的有效性,选取仿真风洞条件下不同高度、不同马赫数的高速飞行器目标图像数据集进行仿真实验,同时选取文献[14]、文献[15]方法进行对比。数据集图像背景为大气背景,包含强噪声、弱目标强度场景;数据集图像中目标为飞行器,可见区域为飞行器与流场摩擦产生的强辐射(飞行器本身不可见)。

同时,为了进一步验证本文算法的性能,对复原图像进行客观质量评价。本文使用峰值信噪比PSNR、结构相似性SSIM、相关系数CORR2、梯度相似性GSM[21]、图像质量指数UQI[22]全面评价图像质量,各指标均为数值越大,图像质量越高,具体结果如图5～图9,表1、2所示。

表1 数据集1复原图像平均质量评价

表2 数据集2复原图像平均质量评价

图5 数据集1样本1图像复原结果

图6 数据集1样本2图像复原结果

图7 数据集1样本3图像复原结果

图8 数据集2样本1图像复原结果

图9 数据集2样本2图像复原结果

图10 复原结果PSNR指标

数据集1为高度30 km、马赫数为5、温度为300 K的长曝光退化图像;数据集2为高度30 km、马赫数为8、温度为500 K的长曝光退化图像。可以看出本文提出的综合图像复原方法,能够良好地保持图像边缘细节,从视觉效果上更清晰,从质量评价指标来看,峰值信噪比PSNR能够达到最优,且在不同噪声强度下,其余评价指标能够平均稳定最优,这说明本文的复原效果较好,且在不同噪声强度、不同飞行条件下都相对更优;文献[33]方法在噪声较强时仍存在偏向某一方向的伪影,文献[34]方法会损失图像的边缘信息,导致SSIM数值较低,抗噪性能也较差。综合来看,本文方法复原结果有更强的目标显著性,在边缘轮廓上更加清晰,整体复原效果更好。

5 结 论

本文针对大气层、临界层背景高速飞行条件下飞行器退化图像,提出了一种面向长曝光图像的图像复原方法,构造多尺度多方向形态学滤波器进行去噪预处理,使用基于图像分块的Radon变换法去除运动模糊伪影,提出基于图像梯度的对称频谱估计法重建气动光学点扩散函数。与其余图像复原算法相比,本文能够自适应确定参数,无需人机交互,同时提升图像细节的分辨能力,复原图像边缘更加清晰,目标图像显著度明显提升,视觉效果更好,可适用于多种类型数据集,满足实时复原的场景。