基于维纳滤波和综合评价因子的遥感图像复原①

王允森，王 勇，左 晨，孟垂哲

(山东航天电子技术研究所，烟台 264670)

0 引言

在遥感图像处理领域，图像的复原和重建是一个非常重要的过程。通常来讲，遥感图像的退化是由各种因素引起的，如大气扰动、拍摄震动、目标辐射、电子信号转换等[1-2]。遥感图像复原的目标即是将退化的遥感图像尽量恢复到原始状态。

传统的图像复原模型一般采用逆滤波方法[3]。由于大部分的图像信息都存在于图像中的边缘，因此，采用滤波方法进行图像复原不仅要去除大部分模糊和噪声，还要保留图像的细节。由于图像的细节和噪声经常在同一区域，这使得图像的去燥平滑和保留细节变为一对相互矛盾的问题。传统的方法难以解决这个问题[4-5]。为了解决这个矛盾，许多改进的逆滤波的方法提出，包括正则化方法、EM算法[6]和GEM算法[7]等。目前，常用的遥感图像半盲复原或盲复原算法一般基于点扩散函数展开，比如Richardson-Lucy方法、维纳滤波方法、最小二乘方法[8-9]等。

近几年，随着处理器计算能力的不断增强，许多学者的研究重点转移到神经网络方法上。这些方法采用基于模型的神经网络实现遥感图像幅复原[10]。Zhou 提出的神经网络复原方法采用Hopfied神经网络来控制图像复原的效果[11]。但是，用Zhou的方法存在网络结构大、网络状态更新慢等问题。为了克服这些缺点，Paik等提出了一种改进的Hopfied神经网络灰度遥感图像复原方法[12]，此后，该算法又被诸多学者再网络模型上进行了改进[13-16]。但是这类算法在计算权重矩阵是都需要耗费大量的空间和时间。后续的深度学习方法大都在此方法的基础上，尽量减少运算时间和运算空间。

另一种广泛应用的图像复原方法是经典维纳滤波算法。标准的维滤波器是一个线性的、空间不变的滤波器。维纳滤波有很多不同的表达方式，信噪比模型，统计模型、参数估计统计模型等[17]。这些方法都可以称之为“维纳滤波”，但是这些方法的原理却不尽相同。采用维纳滤波进行遥感图像复原的关键问题是如何计算原始图像和退化图像的功率谱比，传统的方法一般采用常数代替，这种方法一般会引入新的噪声。

为了克服传统维纳滤波方法的不足，在传统维纳滤波的基础上提出了一种基于改进维纳滤波和综合评价因子的遥感图像复原算法，该算法采用维纳滤波进行退化光学图像的复原，并采用直线运动方法估算图像的退化函数。噪声图像和原始图像的功率谱比采用迭代方法计算，而迭代的依据为综合评价因子。该评价因子集成了三种无参考质量评价方法，包括标准差、对比度和空间频率。为了进一步提高遥感图像质量，采用BM3D方法对维纳滤波复原后的图像进行二次降噪处理。实验表明，提出的图像复原算法能够取得良好的遥感图像复原的效果，并能将运算时间控制在可接受的范围内。

1 遥感图像退化和复原模型

在进行遥感图像复原前，需要根据图像的参数建立的图像退化模型，该模型可以反映出遥感图像退化的原因和机制[18]。如图1所示，f(x,y)是原始图像，g(x,y)是退化的图像，h(x,y)是点退化函数，n(x,y)是高斯噪声函数。

图1 图像退化和复原模型

假设是h(x,y)一个线性空间不变系统，则遥感图像的退化过程可以用式(1)表示

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

(1)

式(1)的傅里叶变换函数可以表示为：

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

(2)

通过式(2)可以算出原始图像的傅里叶表达式：

(3)

假设r(x,y)表示复原后的图像函数，则遥感图像复原的目的则是使下式成立。

R(u,v)=F(u,v)

(4)

2 改进维纳滤波和综合评价因子的遥感图像复原

2.1 维纳滤波

采用维纳滤波方法进行图像复原是建立在系统退化函数和噪声统计的基础上。通常情况下，退化后的图像和噪声函数都是随机函数。而维纳滤波的目标是找到一个最优的滤波器，使得复原后的图像和原图像之间的差值最小。

根据遥感图像复原模型和维纳滤波的复原机制，图像复原问题可以转化为函数的最小化问题，如式(5)所示。

(5)

式中，r(x,y)为复原后的图像，e2表示复原后的图像与原始图像之间的均方误差，E{·}表示数学期望。

当计算均方误差时，频域空间的复原图像函数可以通过式(6)计算。

(6)

式中，H(u,v)是退化函数，H*(u,v)表示H(u,v)的共轭复函数，Sn(u,v)=|N(u,v)|2表示噪声函数的功率谱，Sf(u,v)=|F(u,v)|2表示退化前图像的功率谱。

因此，维纳滤波的转化方程可以表示为：

Hwiener(u,v)=

(7)

除非u=v且H(u,v)=Sn(u,v)=0，否则退化函数不会等于0。维纳滤波方法是建立在最小统计值的基础上。所以从某种意义上说，维纳滤波是一种遥感图像复原的最优算法。

当退化图像的信噪比相对比较高时，Sn(u,v) ≤Sf(u,v)，此时，Sn(u,v)/Sf(u,v)的值将会非常小。这时，维纳滤波变为传统的逆滤波方法，如式(8)所示。

(8)

另一方面，Sn(u,v)≥Sf(u,v)时，维纳滤波器Hwiener(u,v)→0。此时，维纳滤波能够避免逆滤波在图像复原过程中噪声函数过于放大的问题。

在进行维纳滤波之前，需要获得未退化图像的功率值Sf(u,v)和噪声函数的功率值Sf(u,v)。在实际应用时，这两个值一般都是未知的，也难以准确估算。目前常用的方法是令Sn(u,v)/Sf(u,v)等于一个常数矩阵C。这样，维纳滤波的表达式改为：

(9)

式中，C表示未退化图像和噪声图像的功率比。C的值对于图像复原的影响非常大。当C比较小时，意味着信噪比非常小，维纳滤波复原后的图像轮廓会更加清晰，图像更加明亮，但是对噪音的控制非常有限。另一方面，如果C比较大，则复原后他图像能够在一定程度上抑制噪声的放大，但是图像的清晰度往往较差。

因此，采用维纳滤波进行遥感图像的复原需要解决两个关键的问题：退化函数h(x,y)的估计，功率谱参数C的计算。

2.2 退化函数估计

根据遥感图像退化的模型可以得出，如果退化图像的点扩散函数或退化系统的传递函数能够估算出来，则维纳滤波或逆滤波可用于处理遥感图像的退化。这里算法中主要考虑匀速直线运动造成的遥感图像退化。假设退化的方向已知，则退化的图像可以通过旋转使得退化方向转到水平方向上，二维的点扩散函数估算问题可以变成一维[19]。退化系统的点扩散函数满足下式。

(10)

式中，d式退化尺度。这里采用自相关函数来计算退化尺度d。具体步骤如下：

步骤1：在垂直方向对退化图像g(x,y)进行差分运算。

gy(i,j)=g(i,j)-g(i,j-1)

(11)

步骤2：在水平方向对式(11)进行差分运算。

gyx(i,j)=gy(i,j)-gy(i-1,j)

(12)

步骤3：在式(12)的基础上，计算退化图像的每一行的自相关函数。如果第j行有l个像素点，则自相关函数可以表达为：

(13)

式中M是图像的宽度。

步骤4：计算所有自相关函数的均值，并画出自相关函数的曲线。通过寻找曲线的最小值，便可计估算出退化图像的退化尺度。

2.3 功率谱比C估算

对于退化的遥感图像的复原，并没有相对应的清晰图像进行对比。所以，这里采用无参考图像质量评价方法对遥感图像进行复原参考，这里采用了3种常见的图像质量评价方法：灰度标准差、对比度、空间频率[20-22]。

(1)标准差(Standard Deviation,STD)

标准差STD表示了图像像素灰度值相对于均值的离散程度。如果标准差越大，表明图像中灰度级越分散，图像质量也就越好。标准差是最简单的评估图像质量的方法，用来评估图像灰度的变化，其计算公式为：

(14)

式中，M,N分别表示遥感图像的宽度和高度，u表示图像像素的均值。

(15)

(2)对比度(Contrast Ratio,CR)

图像的对比度代表了图像黑与白的比值，也就是遥感图像从黑到白的渐变层次，比值越大，从黑到白的渐变层次越多，图像的色彩越丰富。采用四邻域方法计算图像的对比度，首先对原始图像的四个边界进行复制扩展，假设图像由f变为fs，而图像的大小由M×N变为(M+2)×(N+2)，对比度的计算公式如下：

(16)

其中，δ(i,j)，表示图像fs相邻像素间的灰度差，计算方法如下：

δ(i,j)=[fs(i,j)-fs(i-1,j)]2+[fs(i,j)-fs(i,j-1)]2

+[fs(i,j)-fs(i+1,j)]2+[fs(i,j)-fs(i,j+1)]2

(17)

(3)空间频率(Spatial Frequency,SF)

空间频率SF指标反映了融合图像在空域的总体活跃程度，空间频率值越大，对应融合结果越好，具体定义为：

(18)

其中，RF和CF分别为图像的行频率和列频率。对于大小为M×N的待评价图像f(x,y)，其行频率和列频率计算如下：

(19)

(20)

为了对遥感图像进行综合评价，提出的综合评价因子Qc通过STD、CR和SF计算。

(21)

式中，w1,w2和w3分别表示三种评价方法的权重参数。STDr,CRr和SFr是复原后图像的3个评价参数，STDg,CRg和SFg是复原前图像的3个评价参数。

获得最佳功率谱比C意味着提出的综合评价指标Qc必须越大越好。

这里采用一种自迭代维纳滤波方法计算功率谱比C的最优值。步骤如下

步骤1：设置功率谱比C0、迭代步长ΔC和迭代次数Ct的初值，令C=C0；

步骤2：通过功率谱比的值，对遥感图像进行维纳滤波，并计算图像的综合评价因子Qc。

步骤3：更新功率谱比C的值，C=C0+ΔC.迭代次数t=t+1.

步骤4：循环执行步骤2和步骤3，直到t=Ct.经过迭代，获得了一系列图像综合评价因子的值。

步骤5：通过Qc的最大值，找到最优的功率谱比Copt。如果Copt不是C的第一个值或最后一个值，则Copt为真正最优的功率谱比，迭代结束。

步骤6：如果最大的综合评价因子Qc对应C的初始值，则ΔC=ΔC/10，并循环执行步骤2～步骤5； 另一方面，如果Qc对应C的最后一个值，则ΔC=ΔC×10，并重复步骤2～步骤5直到迭代结束。

执行步骤6的主要原因是防止不适当的功率谱比初值C造成局部最小的问题。

3 实验结果及复原效果分析

3.1 仿真验证

在本章中，采用仿真实验来验证提出的改进维纳滤波算法和综合评价因子，测试图像如图2(a)所示，大小为512像素×512像素。实验过程首先对该图像进行匀速直线运动的PSF退化处理和高斯噪声污染，如图2(b)所示，退化尺度设为5，迭代参数C0=0,ΔC=0.001,Ct=100，为了确定三个评价因子的权重系数，分别采用三种评价方法作为单独对测试图像进行复原，并统计综合评价因子的大小。

(a)原始图像 (b)退化图像

采用标准差、对比度、空间频率作为单独评价因子，分别对退化图像进行复原，复原效果如图3(a)～(c)所示，对应的综合评价因子和功率谱比的值如表1所示。可以看出，采用标准差进行图像复原的参考因子时，复原图像中的前景与背景区分比原图和退化图像更加明显，但图像中的某些细节被覆盖掉。图像的对比度代表了图像黑与白的比值，也就是图像从黑到白的渐变层次，然而退化图像中包含了许多噪声信息，所以只采用对比度作为单独评价因子时，复原结果的渐变层次反映在了噪声上，处理结果不理想。而空间频率关注图像上下文纹理的丰富程度，复原图像比较接近原始图像，但是背景和前景的对比要弱于标准差。

表1 三种评价因子统计比较

遥感图像的处理一般是为了使得图像的背景和目标对比更加明显，但需尽可能的保留图像纹理信息。根据三种评价因子的复原结果和统计值，在进行遥感图像复原时，综合评价因子中，应当使得标准差所占比重最大，对比度所占比重最小。

为了进一步说明改进维纳滤波算法的优点，将传统维纳滤波方法和逆滤波方法同改进的维纳滤波方法进行了对比。根据以上分析，STD、CR、和SF三个因子的比重分别设置为0.5、0.2、0.3，测试对比结果如图4所示。

(a)逆滤波 (b)传统维纳滤波 (c)改进的维纳滤波

对比图4中三种方法的复原结果，从人眼的观察角度而言，改进的维纳滤波方法的复原效果比其他两种方法都要好。文献[9]指出，如果要改图像未被噪声污染，则逆滤波方法是一种比较好的复原方法。但是，在存在噪声污染的前提下，逆滤波的效果难以令人满意，传统的维纳滤波方法同样如此。

传统维纳滤波执行过程中，功率谱比CSNR通过图像的信噪比计算，即CSNR=N(u,v)/F(u,v)。

图4(b)显示了传统维纳滤波的复原效果。可以看出，由于新的噪声的引入，复原后的图像质量甚至不高于污染后的图像图4(b)。表2统计了原始图像、退化图像和三中复原方法的评价参数计算结果和运算时间。可以看出，提出的改进维纳滤波算法在存在高斯噪声和直线退化的影响下，复原效果相对有效，用该方法计算的峰值信噪比PSNR也相对较大。此外，由于改进的维纳滤波需要多次复原图像，并迭代计算综合评价因子大小，因此运算时间较逆滤波和传统维纳滤波有所提高，对于512×512的测试图像，运算时间为141ms，是逆滤波的2.5倍，传统维纳滤波的1.7倍，并不是数量级的提升。但综合评价因子大小却是另外两种方法的8.2倍、19.2倍，图像质量有了显著提高。

表2 三种处理方法比较

3.2 遥感图像复原结果

为了进一步说明改进维纳滤波的优势，选择4张大小为2500像素×2500像素的遥感图像进行算法验证，如图5所示。这些遥感图像都经过直线点扩散函数的退化和高斯噪声的污染，但是退化的尺度未知。根据遥感图像的处理经验，初始的参数设置如下：C0=0.0025, ΔC=0.001,Ct=100，三种评价参数的权重分别设为0.5，0.2和0.3。此外，遥感图像往往存在非平稳的随机过程的噪声，采用成熟的BM3D方法对维纳滤波复原后的图像进行二次降噪处理。

(a1)原始遥感图像 (b1)维纳滤波复原后的遥感图像 (c1)BM3D方法去噪后的图像

(a2)原始遥感图像 (b2)维纳滤波复原后的遥感图像 (c2)BM3D方法去噪后的图像

(a3)原始遥感图像 (b3)维纳滤波复原后的遥感图像 (c3)BM3D方法去噪后的图像

(a4)原始遥感图像 (b4)维纳滤波复原后的遥感图像 (c4)BM3D方法去噪后的图像

从图5的原始图像和复原图像的比较中可以看出，复原图像的轮廓和图像细节更加清晰。图像中的地面特征更加明显和明亮。表3为4幅遥感图像三种评价因子的大小。可以看出，采用基于综合评价因子的改进维纳滤波进行图像复原，图像的质量有了一定程度的提高。而经过BM3D算法后，遥感图像的噪音进一步减少，可作为目标检测、识别的输入对象。

表3 处理前后质量参数对比

为了进一步测试算法在星载平台的处理性能，选择嵌入式图像处理运行平台测试算法运行时间。Nvidia Tegra X2是目前广泛应用的嵌入式GPU平台，功耗仅为7.5瓦，浮点运算能力能够达到1.3Tflops，目前在工厂机器人、无人机和智能摄像头等领域应用广泛，并已在某些实验卫星成功应用，主要用于遥感图像的快速预处理和目标检测识别等。因此选择Nvidia Tegra X2作为遥感图像处理平台。经统计，算法在Nvidia Tegra X2平台复原100副相同大小的遥感图像的平均时间为4.54秒，对于2500×2500大小的图像是可以接受的。

4 结束语

通过分析光学遥感图像退化的原因和维纳滤波的原理，提出了一种新的基于改进维纳滤波和综合评价因子的遥感图像半盲复原算法，与传统的维纳滤波方法相比，提出的方法具有以下优势：(1)通过基于综合评价因子的迭代方法，可以获得更加精确的功率谱值，这对遥感图像的复原相当重要；(2)提出的半盲复原算法集成了对比度、标准差和空间频率三种无参考图像质量评价因子，这三种因子在一定程度上代表了原始未退化图像的属性。实现表明提出的算法能够获得良好的遥感图像复原效果，提出的算法能够对遥感图像进行复原处理。