基于薄云厚度分布评估的遥感影像高保真薄云去除方法.docx

该【基于薄云厚度分布评估的遥感影像高保真薄云去除方法 】是由【科技星球】上传分享，文档一共【16】页，该文档可以免费在线阅读，需要了解更多关于【基于薄云厚度分布评估的遥感影像高保真薄云去除方法 】的内容，可以使用淘豆网的站内搜索功能，选择自己适合的文档，以下文字是截取该文章内的部分文字，如需要获得完整电子版，请下载此文档到您的设备，方便您编辑和打印。基于薄云厚度分布评估的遥感影像高保真薄云去除方法
 
 
汪月云,黄 微,王 睿
(上海大学 通信与信息工程学院，上海 200444)(*通信作者电子邮箱******@)
0 引言
光学遥感影像往往受到薄云污染的影响，造成影像地物细节对比度下降，这给影像的应用[1-3]带来了很多困难。因此，如何有效去除降质影像的薄云影响从而提高影像的利用率是遥感影像应用的一项重要研究。
目前，光学遥感影像薄云去除主要有频域法和空域法两类。频域法去云主要根据薄云信息处于影像低频区域的原理，所以可以设计合适的高通滤波器以滤除低频成分从而有效去云。比如:阎庆等[4]结合小波变换和同态滤波、冯春等[5]利用同态滤波都达到了有效去云效果。
空域法去云将薄云污染分为乘性污染和加性污染两种模型。乘性污染经典的去云方法有暗通道先验法[6-7]和薄云优化变换(Haze Optimized Transformation, HOT)法[8-9]。暗通道先验法通过对清晰影像库统计得到暗通道分布规律，即清晰影像至少在某一个颜色通道中有灰度值接近于0的像素，利用该先验规律得出云模型的大气光和透射率从而进行薄云去除。Pan等[10]研究发现无云遥感影像的暗通道数值很低，但不接近于0，通过改进遥感影像暗通道分布从而有效去除薄云。HOT法是利用遥感影像红光波段和蓝光波段相关性高的特点寻找合适的晴空线，用像元到晴空线的距离表征该点云量进而有效去云。He等[11]提出了AHOT(Advanced HOT)算法，先利用空间信息剔除易混淆地物,提高云检测精度，再考虑大气散射乘性效应，提出基于AHOT的虚拟云点(Virtual Cloud Point, VCP)技术扩大由薄云污染降低的数字量化值(Digital Number, DN)方差从而有效去云。Makarau等[12]利用加性薄云污染模型，在传统的暗目标减法(Dark Object Subtraction, DOS)基础上[13]，通过局部搜索整个影像的暗目标构建影像的薄云厚度分布(Haze Thickness Map, HTM)，最后用原降质影像减去薄云厚度分布得到无云影像。这些去云算法大多关注去云后影像细节的锐化、地物间对比度的增强，如何在去云的同时保证数据保真是一项具有挑战性的问题。
近年来，已有学者提出了一些考虑保真的薄云去除算法[14-16]。Shen等[14]基于同态滤波法半自动化地确定不同波段的截止频率，对薄云像素和无云像素分别处理从而获得高保真的无云影像。Singh等[15]基于暗通道先验法提出改进的无云影像恢复算法，通过改进联合三边滤波器以提高大气光的精细估计效果，并重新定义透射率从而降低无云影像的颜色失真程度。Liu等[16]基于加
性薄云污染模型，将原始薄云厚度分布先经过纹理去除，再进行大面积地物反射抑制，然后用导向滤波器进行平滑，最后用线性回归求取不同波段的薄云厚度从而得到无云影像。文献[16]算法的去云结果整体具有保真效果，但是依旧存在无云区域过度校正，以及影像中蓝色地物色彩严重失真的问题。
因此，本文在文献[16]算法的基础上提出了一种改进的薄云去除方法，即基于薄云厚度分布(HTM)评估的高保真薄云去除方法。该方法首先利用基于暗通道先验薄云去除算法求得影像的透射率分布图[10]，根据透射率分布图中像素数值接近于1来自动确定无云区域的范围。然后,在HTM中，求取无云区域的平均值,将整个HTM减去无云区域平均值，使得HTM满足无云区域的薄云厚度接近于零，以确保无云区域的保真度。针对蓝色地物区域云量被过估计，造成去云结果中色彩失真的问题，本文利用蓝色地物RGB波段的光谱特性，结合阈值法检测出蓝色地物，用蓝色地物附近的HTM平均值作为蓝色地物区域的云量估计。最后,用降质影像减去不同波段优化的HTM得到最终去云结果。
通过多幅影像的实验发现，与文献[14]算法去云结果、文献[15]算法去云结果，以及文献[16]算法去云结果相比，本文方法的结果影像去云效果良好，同时具有数据高保真的优点。
1 加性模型薄云去除
云成像模型
薄云污染可以认为是加性污染[12，16],并且薄云污染对不同波段的影响不同[16-17],所以云成像模型可以描述为：
Ic(x)=Jc(x)+Hc(x);c∈{R,G,B}
(1)
其中：x代表像素的位置；I代表降质影像；J代表无云影像；H代表薄云厚度分布；c是{R,G,B}中的任一波段。
薄云估计
在Makarau等[12]的研究中，薄云厚度分布(HTM)估计如下：
(2)
其中：Ib(x)、Ig(x)分别是降质影像的蓝光波段和绿光波段；W(x)是以x为中心、大小为w×w像素的窗口。本文实验中该窗口大小取3×3像素。
Liu等[16]发现由式(2)得到的Hraw(x)存在将地物反射光作为薄云组成成分，造成降质影像被过度矫正的问题。为了减少地物反射光的影响，文献[16]对Hraw(x)进一步改进：首先为了去除Hraw(x)中的地物纹理，文献[16]采用了可以去除无用纹理和细节的全变差(Total Variation, TV)正则化图像修复模型[18-19]得到平滑的薄云厚度分布Hs(x)；接着,为解决大面积高亮地物的薄云厚度过估计问题，文献[16]用Cannny算子检测出高亮地物边缘，用局部窗口对Hs(x)中高亮地物区域作数值抑制得到Hd(x)；然后,应用导向滤波器[20]去除Hd(x)中的块状效应得到Hg(x)；最后,求取不同波段的薄云厚度分布。
Himpr(x)=Hg(x)-minHg
(3)
Hc(x)=kcHimpr(x)；c∈{R,G,B}
(4)
为了去除气溶胶的影响，式(3)用Hg减去其最小值minHg得到最终的薄云厚度估计Himpr(x)。式(4)中：kc是不同波段的云估计强度，它是在Hg低于平均值的像素集与降质影像I的不同波段中低于平均值的像素集的交集中，以Hg的像素为自变量，以影像的不同波段的像素为因变量，通过线性回归得到；Hc(x)是不同波段的薄云厚度估计。得到Hc(x)后，通过式(1)即可得到无云影像J。
图1是根据文献[16]算法得到的实验结果。其中:图1(a)是大小为300×400像素的GF-2降质影像;图1(b)是最终的云厚度分布估计Himpr(x);图1(c)是去云结果。由图1可以看出，Himpr(x)中尽管云厚度分布范围估计比较准确，但是无云区域的云厚度值并不为0或接近于0，同时无云区域的蓝色和白色物体估计的云厚度分布值也很高，这就造成了去云影像的无云区域存在过度矫正和蓝色地物区域严重失真的问题。
图1 文献[16] Haze estimation and haze removal results of algorithm in literature [16]
2 薄云去除改进
为解决图1去云结果中反映出来的无云区域过度矫正和蓝色地物失真严重问题，基于加性薄云污染模型，本文在文献[16]得到的最终薄云厚度分布Himpr(x)基础上进行改进，在获得良好去云结果的同时确保影像的高保真。本文方法的流程如图2所示。
图2 本文方法流程Fig. 2 Flow chart of the proposed method
无云区域保真
首先针对无云区域薄云厚度被过估计，造成无云区域数据失真的问题，本文对Himpr(x)作了如下改进：
Hnew(x)=max(Himpr(x)-K,0)
(5)
其中：
(6)
式中：Dclean是Himpr(x)中的无云区域;M是Dclean中的像素个数;K是无云区域的平均值。理论上，Hnew(x)的无云区域DN值接近或等于0，有云区域DN值大于0，为了避免由于Himpr(x)-K导致负数的出现，这里使用了最大值算子max。
无云区域Dclean的确定影响着最终无云影像无云区域的保真度。为了有效确定无云区域Dclean的范围，本文利用基于暗通道先验薄云去除算法[6-7]中的透射率来估计该区域。其中透射率为：
t(x)∈(0,1]
球形工频电场检测传感器采用球形电容检测感应电压的方法来检测电场,一维球形传感器的结构如图1所示,将一个空心金属球壳切割成两部分,通过绝缘物质将两个半球连接在一起,并在两个半球中间并联一个测量电容,两个半球分别构成了传感器的两个电极[6-7]。
(7)
蓝色地物云量估计
对于降质影像中的蓝色地物，由于像素在蓝光波段的亮度值远大于绿光波段，因此根据式(1)得到的薄云厚度分布Hraw(x)中容易造成蓝色地物薄云厚度被过估计，最终导致去云结果中蓝色地物色彩失真。
通过大量实验数据发现，在优化后的薄云厚度分布Hnew(x)中，白色地物和蓝色地物区域的DN值都很高。因此，本文先利用阈值选定白色区域和蓝色区域，再根据原始影像中蓝色地物蓝光波段像素值远大于红光波段像素值的特点[21]，剔除白色地物区域最终得到蓝色地物的像素集：
P1={x|Hnew(x)≥T}
(8)
P2={x|Ib(x)-Ir(x)>0}
(9)
Pf=P1∩P2
(10)
其中：x代表像素坐标位置；P1是包含蓝色地物和白色地物的像素集；P2是包含偏蓝地物的像素集；Pf是蓝色地物像素集。通过阈值T可以选定白色地物区域和蓝色地物区域，阈值T的选定不能过大，否则影像中的蓝色地物区域不能全部选中，阈值T也不能过小，否则容易包含非蓝色地物。本文中阈值T取Hnew(x)中最大像素值的70%，即T=max(Hnew(x))*。Ib(x)、Ir(x)分别是原始影像的蓝光波段和红光波段。由于薄云的分布一般是连续和平缓的，因此本文将蓝色地物区域的薄云厚度用周围区域薄云厚度的均值代替，其表达式为：
(11)
其中:Hblue是蓝色地物薄云厚度;Hnew(x)是优化的薄云厚度分布;Dnear是在Hnew(x)中蓝色地物附近人工选取出的薄云区域;n是Dnear中像素点个数。
经过以上的优化，可得到完整的薄云厚度分布Htotal(x)：
(12)
为了去除Htotal(x)中的块状效应，本文应用导向滤波器[20]将其平滑得到最终的薄云厚度分布Hfinal(x)。
本文方法的每一个优化步骤结果如图3所示。其中：图3(b)是根据透射率得到的无云区域二值图像，白色表示无云区域。图3(c)是优化后的薄云厚度分布Hnew(x)。可以看出，Hnew(x)在无云区域的DN值接近或等于0，这能保证恢复的无云影像在无云区域的高保真度。图3(d)是蓝色地物的检测结果, 可以发现，图中的蓝色地物位置检测比较准确。图3(e)是最终的薄云厚度分布Hfinal(x)。与图3(c)相比，图3(e)的Hfinal(x)在蓝色地物区域的DN值相比Hnew(x)小很多，可以避免蓝色地物云量被过估计。
图3 本文方法各步骤结果Fig. 3 Results of each step of the proposed method
高保真无云影像
通过逐步优化薄云分布，单独估计蓝色地物云量，最终得到的无云影像J表示如下：
Jc(x)=Ic(x)-kcHfinal(x)
(13)
其中：I是降质影像；kc是根据文献[16]算法中通过线性回归得到的不同波段云估计强度；Hfinal(x)是最终优化的薄云厚度分布。图3(f)是本文改进方法得到的无云影像，对比图3(a)的降质影像可以发现，图3(f)的无云影像不但具