2. 安徽大学 电子信息工程学院，安徽 合肥 230039
2.School of Electronics and Information Engineering, Anhui University, Hefei 230039, China
1 引言
遥感影像在获取时,云层会遮挡地物,也会在影像中产生阴影。去除影像中的云层及云阴影,再现清晰地物,以便于遥感影像的进一步分析处理和应用,具有十分重要的意义。
小波变换去云方法是目前使用较多的厚云去除方法[1,2]。该方法利用小波变换的多尺度、多分辨率特性,依据一定的灰度值修正范围进行影像融合达到去云的目的。其缺点是两幅影像的厚云区域要求互不重叠,多数情况下不能满足实际需要。
影像融合去云方法[3]将同一地区的多源或多时相遥感影像进行融合,达到去云目的。但这需要多幅影像的云区不重叠,而且需要解决可能存在的辐射差异问题。
图像修复法[4,5]利用云层和云阴影所在区域、与清晰地物的过渡区域以及清晰地物区域的统计学特性,对图像中的云层和阴影灰度值进行修正实现去云。这种方法的优点在于只需要单幅图像即可实现去云,缺点是只能针对云层和阴影区域范围较小的情况。
此外还有多光谱图像法[6],直方图匹配法,同态滤波法,基于学习的方法[7]等。这些厚云去除方法中,或需要图像的特殊光谱信息,或要求图像纹理结构强,且厚云的去除不彻底,影响了遥感影像的正确判读。
本文提出一种基于支持向量机的遥感影像厚云及云阴影去除方法,利用支持向量机优异的学习性能,并结合太阳方位角和高度角信息,有效地检测厚云及其阴影,通过图像镶嵌和高低频补偿的方法去除影像中的厚云及云阴影,与小波变换等方法相比,能更好地恢复厚云区域的地物信息。
2 支持向量值轮廓波变换 2.1 支持向量值滤波器文献[8]利用最小方差支持向量回归模型提出了支持向量值滤波器。最小方差支持向量回归模型可表示为
式中,K={Ki,j=K(xi,xj)}Ni,j=1;e=[1 … 1]T;Y=[Y1 … YN]T;α=[α1 … αN]T;γ是均衡常数。
定义
则式(1)的解为
将模型式(1)应用于图像处理中,其中图像块的像素坐标值作为训练样本的输入值,即式(1)中核矩阵K的输入,图像块的像素值作为训练样本的输出值,即式(1)中的Y值,拉格朗日乘子α的值可以看做处理后的分解系数。
对于任意图像窗口,像素坐标具有下述形式
通过减掉中心坐标(x0,y0),所有这样的点集可以变换成统一形式
对于任意图像窗口,核矩阵K以及式(3)中的矩阵Q就有相同的取值。将矩阵Q的中心行向量排列成方阵,就获得了支持向量值滤波器[8]。 2.2 支持向量值轮廓波变换
设第j尺度下的低频图像为Pj,高频细节图像为Sj,j=1,2,…,r,第j尺度下的支持向量值滤波器为SVj,支持向量值变换的分解过程可以表示为
方向滤波器组对高频细节图像Sj进行方向滤波过程可以表示为[9]
式中,S(lj)j,k是方向滤波器组在第j层的第k个分量。
支持向量值轮廓波变换的重构过程为[9]
3 厚云及云阴影去除方法
首先进行厚云及云阴影检测,利用多时相遥感影像,通过影像镶嵌方法去除影像中的大片云层及云阴影,再对小范围的云层及阴影重叠区域通过统计学补偿的方法进行修复,完成厚云及云阴影去除[10]。
3.1 影像预处理影像预处理主要包括地理位置和亮度的配准。地理位置配准是通过选择一定的地面控制点,利用多项式模型方法配准影像。亮度的配准是使用线性扩展方法,利用下式对待纠正影像进行调整,得到与基准影像相同灰度范围的影像[11]
式中,[a,b]和[c,d]分别是基准影像和待纠正影像的灰度范围。 3.2 云层及云阴影检测
从原始影像中提取部分像素点作为有云区和无云区的样本,利用支持向量机进行分类[12],生成一幅二值图,图中的黑色区域为无云区域,白色区域为云层区域。
云阴影检测的方法如图 1所示,利用太阳方位角和太阳高度角,并结合云层高度,可以确定云阴影相对于云层的方位以及云阴影边缘相对于云层边缘的实际位移,从而预测阴影区域,得到阴影的二值图[13]。
|
| 图 1 云阴影检测示意图 Fig. 1 The diagram of cloud shadow detection |
设影像A与影像B是需要镶嵌的两幅含云遥感影像,假设A是基准影像,B是参考校正影像,把B中的一部分ΩB镶嵌到A中,如图 2所示。
|
| 图 2 影像镶嵌示意图 Fig. 2 The diagram of image mosaic |
令
式中,Sj,kA(x,y)为影像A的第j层第k个方向像素坐标为(x,y)的子带系数;PjA(x,y)是影像A的第j层像素坐标为(x,y)的低频系数;CjΩ(x,y)定义为
对调整后的高低频系数进行式(9)的重构,即可得到镶嵌影像[14]。
3.3.2 云层初去除为了解决两幅遥感影像云层互相重叠的问题,利用多源影像进行多幅影像镶嵌,下面是云层初去除的算法步骤:
(1)进行基于支持向量机的云层检测,对每幅影像生成云层二值图,其中有云区域值为1,无云区域值为0。
(2)利用云层二值图,构造图像的一级选择矩阵。以3幅图像为例,分别记为图像1、图像2和图像3,式(14)说明了图像一级选择矩阵的构造方法。
式中,1-cloud、2-cloud、3-cloud分别为图像1、图像2和图像3的云层二值图;array1、array2、array3、array4为图像一级选择矩阵。按式(14)赋值,其余元素为零。
(3)利用支持向量值轮廓波变换对源影像进行分解,将分解后的图像2和图像3分别乘以相应的一级选择矩阵array2、array3,将非零结果替换图像1中相应区域,实现云层的初步去除。
利用式(14)的一级选择矩阵,图像1中的非云区将保持不变,若某区域是图像1的云区且是图像2的非云区,则用图像2中该区域镶嵌进图像1中;若其是图像1和图像2的云区且是图像3的非云区,则用图像3中该区域镶嵌进图像1中;若其位于3幅图像的云层重叠区,则留给后节进行补偿处理。记云层初去除的图像为cloudremoval_1。
通过上述云层初去除方法,实现了云层初去除的"多层镶嵌",如图 3所示。
|
| 图 3 云层初去除的“多层镶嵌” Fig. 3 “Multilayer mosaic” of preliminary cloud removal |
这种"多层镶嵌"使得逐步镶嵌的非云区域不相关。为了去除"低阶"影像中的云层,只能通过逐步镶嵌"高阶"影像中相同位置的非云区域,这里所说的"低阶"和"高阶"是根据嵌入的先后次序来分的,即需要去云的基准影像处于最"低阶",最后用来镶嵌的参考影像处于最"高阶"。
3.3.3 云阴影初去除完成云层初去除后,需要进行云阴影初去除。云阴影初去除也是通过设计二值图来进行,过程如下:
(1)利用太阳角度信息进行云阴影检测,生成云阴影二值图,其中阴影区域值为1,无阴影区域值为0。计算cloudremoval_1所含的图像1云阴影部分,记为cr-1。结合云检测二值图,可以得出3幅影像地物二值图,记为1-land,2-land,3-land,其中清晰地物区域值为1,其余区域值为0。
(2)计算图像的二级选择矩阵。式(15)是以3幅图像为例,说明图像二级选择矩阵的构造方法。
这里Array1、Array2、Array3为图像二级选择矩阵,按式(15)赋值,其余元素为零。
(3)将图像2和图像3的分解系数分别乘以二级选择矩阵Array1、Array2,将非零结果替换图像cloudremoval_1中相应的分解系数,实现云阴影初去除。
利用式(15)的二级选择矩阵,对cloudremoval_1中所含的"低阶"的云阴影,利用其"高阶"的清晰地物来镶嵌。若找不到镶嵌的"高阶"的清晰地物,则留待云阴影补偿环节处理。记经过云阴影初去除的图像为cloudremoval_2。
3.4 云层及云阴影的补偿对影像镶嵌未能去除的云层及云阴影,可以通过统计学补偿的方法进行去除。本文对经过云层和云阴影初去除的图像cloudremoval_2的高低频系数分别进行补偿。由于云阴影位于非云区,因而需要先进行云阴影补偿,再进行云层补偿。
3.4.1 云阴影补偿由于图像的低频系数具有慢变化特性,某一区域的低频信息与其周边一定范围内的低频信息具有统计相似性,可以利用此特征进行云阴影的低频系数补偿。
低频补偿按照式(16)的策略进行
式中,P1j(x,y)和LC1j(x,y)分别表示补偿前和补偿后的云阴影的低频系数;mc和σc分别是云阴影周边清晰区域低频系数的均值和方差;ms和σs分别是云阴影区域低频系数的均值和方差。
高频补偿按照式(17)的策略进行
式中,S1j,k(x,y)和HC1j,k(x,y)是补偿前后的子带高频系数;Ωcloud是云区域;Ωclear是清晰地物区域;Ωshadow是云阴影区域;Ωtrans是云阴影和清晰地物之间的过渡区域;a1≥1是高频增强系数[2]。
经过上述高低频补偿并重构后的影像,云阴影已经基本去除,记云阴影补偿后的影像为cloudremoval_3。
3.4.2 云层补偿经过云层初去除之后,只剩下几幅影像的云层重叠区域。对于该云层重叠区域的高低频系数进行补偿,处理方法与阴影补偿相类似。
低频补偿按照式(18)进行
式中,LC1j(x,y)和LCj(x,y)分别表示补偿前和补偿后的云层重叠区域的低频系数;mc、σc是周边清晰区域低频系数的均值和方差;ms、σs是云层重叠区域低频系数的均值和方差。
高频系数补偿如式(19)所示
式中,HC1j,k(x,y)和HCj,k(x,y)分别表示补偿前和补偿后的云层重叠区域的高频子带系数;a2≤1是高频抑制系数。
经过系数补偿后,利用支持向量值轮廓波变换重构高低频系数,得到了经过云层补偿的影像cloudremoval_4。对该影像进行中值滤波平滑,即得到最终的去云结果影像cloudremoval。
4 试验结果分析通过试验分析本文提出的厚云及云阴影去除方法的性能,并与小波变换去云方法、图像融合去云方法和图像修复去云方法进行对比分析。其中小波变换去云方法用的是文献[2]中的方法,图像融合去云方法用的是文献[3]方法,图像修复去云方法首先进行云层检测,再用参考影像的无云像元替换基准影像的有云像元,最后用文献[5]的基于Bandelet的图像修复技术对剩下的小块云区和阴影进行修复。试验中使用的遥感影像为Landsat5 TM拍摄的美国苏瀑布城的3幅多时相遥感影像,见图 4,其中(a)是需要去除厚云及云阴影的基准影像,获取时间是2006年9月9日,太阳方位角是136.25°,太阳高度角是56.11°,(b)和(c)是用来镶嵌的参考影像,影像(b)的获取时间是2006年10月27日,太阳方位角和太阳高度角分别是154.40°、41.88°,影像(c)的获取时间是2006年9月25日,太阳方位角和太阳高度角分别是144.05°、51.74°。
|
| 图 4 美国苏瀑布城TM影像 Fig. 4 TM images of Sioux Falls City of American |
以图 4(a)为基准,将图 4中的参考影像(b)和(c)按本文方法进行预处理,结果见图 5。根据云顶亮温和纹理信息,可以检测出高云、低云等,由于本次试验中所使用的含云遥感影像没有云顶亮温信息,从纹理信息可以判断图中的云层属于低云,低云的高度一般在2000m左右。利用支持向量机检测影像中的云层,并利用太阳角度信息,判定云阴影区域,得到云层和云阴影的二值图,接着使用支持向量值轮廓波变换,通过二值图构造选择矩阵进行影像镶嵌的方法先完成厚云和阴影的初步去除,云层及云阴影初去除结果见图 6,云阴影及云层补偿结果见图 7,小波变换方法、图像融合方法、图像修复方法和本文方法去除厚云及云阴影结果见图 8。可以看出不同的厚云及云阴影去除方法中,本文方法的去云效果比较好。经过一步一步的处理,云层和云阴影逐渐被处理掉,得到去云后的清晰影像,影像中的河流较好地显现出来。小波变换方法由于抑制了图像的低频特性,破坏了图像的概貌信息,引起去云图像的色彩失真。图像融合方法也能有效去除云层和云阴影,但从图 8(b)可以看出,融合结果图仍然带有部分源影像的云层信息,这些云层信息甚至遮挡了图像中央的部分河流。图像修复方法虽然能够大幅消除阴影带来的不良效应,且图像的边缘信息相对小波方法有所提高,色彩保持方面做得也很好,但是,由于依据统计学的原理缺乏必要的先验知识,使得图像中的河流信息未能完全恢复,出现了"断流"现象。本文方法有效去除了云层及云阴影,有效恢复厚云区域的地物信息,形成的无云图像细节清晰,图像光滑。
|
| 图 5 参考影像的预处理结果 Fig. 5 Preprocessing results of reference images |
|
| 图 6 云层及云阴影初去除结果 Fig. 6 Preliminary results of removal of cloud and cloud shadow |
|
| 图 7 云阴影及云层补偿结果 Fig. 7 Compensating results of cloud and cloud shadow |
|
| 图 8 不同方法去除厚云及云阴影结果 Fig. 8 Results of removal thick cloud and cloud shadow for different methods |
均值、影像信息熵、相关系数、偏差指数以及基于视觉感知的影像质量综合评价方法,是本文采用的定量评价标准。
影像去除云层及云阴影后,如果均值比较适中,表明云层及云阴影去除越好,均值偏大或偏小都不是最佳的效果。
影像的熵值是影像信息丰富程度的一个重要评价指标,较大的熵值说明图像的细节表现力较强。
相关系数反映去云图像与源图像的相关程度,偏差指数反映去云图像与源图像的差异程度。
上述评价指标必须综合考虑,才能反映去云效果的好坏。各种方法去除厚云及云阴影的定量评价结果见表 1。
| 原始 影像 | 小波 变换 | 图像 融合 | 图像 修复 | 本文 算法 |
| 均值 | 102.99 | 70.471 | 71.215 | 72.432 | 72.771 |
| 影像信息熵 | 6.844 0 | 5.271 1 | 5.177 6 | 6.056 4 | 5.757 4 |
| 相关系数 | / | 0.080 1 | 0.152 8 | 0.248 3 | 0.217 1 |
| 偏差指数 | / | 0.085 1 | 0.089 6 | 0.082 0 | 0.081 6 |
从表 1可以看出,小波变换方法、图像融合方法和图像修复方法去云后图像的均值要低于本文算法,这是小波变换方法、图像融合方法和图像修复方法对阴影的去除不够完全导致的,图 8(a)、(b)、(c)中间的波浪形痕迹就是没有去掉的阴影。本文算法的影像信息熵和相关系数要高于小波变换方法和图像融合法方法,低于图像修复方法,这是由于小波变换去云方法破坏了图像的低频信息,引起了去云图像的色彩失真,图像融合方法带有部分低频的云层信息,因而信息熵最低,而图像修复方法由于未能完全去除图像中的波浪形阴影,因而从源影像中获取的信息更多。从偏差指数指标来看,本文算法的去云效果也要优于其他方法。
本文还采用一种更为合适的评价标准,即基于视觉感知的镶嵌影像质量评价方法来验证算法的有效性。视觉感知(visual perception,VP)是基于人眼视觉特征的评价标准。可用空间频率活动性(spatial frequency activity,SFA)与结构对比度(structural contrast degree,SCD)指标来评价[15]。
空间频率活动性SFA的计算公式为
式中,F (i,j)是去云影像F在点(i,j)处的灰度值。
结构对比度SCD的计算公式为
式中,s1n、t1n分别是基准影像中待镶嵌区域sn和镶嵌结果影像中镶嵌区域tn像素灰度值的均值;k和l分别表示影像镶嵌前后灰度均值发生变化的像素块的个数。
综合空间频率活动性和结构对比度这两个方面,得到基于视觉感知的综合评价计算公式为[14]
式中,Sx和Sy分别是归一化后的感知因子SFA和SCD。
去云影像的空间频率活动性越小,表明镶嵌边界的拼接缝越光滑,云层及云阴影去除的视觉效果就越好。在结构对比度计算公式中,如果是厚云区域或者云阴影区域,则按前2种情况进行计算,而对于未经处理的非云区域按第3种情况计算。去云影像的结构对比度小,表明影像整体和谐,去云和云阴影效果好,反之,则对视觉刺激大,影像整体和谐性差,去云及云阴影效果不好。综合空间频率活动性和结构对比度两方面指标的视觉感知评价标准,从人类视觉特性上直观地反映了云层及其阴影去除的优劣。基于视觉感知的评价结果见表 2。
| 小波变换 | 图像融合 | 图像修复 | 本文算法 | ||
| VP | SFA | 0.147 9 | 0.079 7 | 0.088 9 | 0.078 1 |
| SCD | 0.687 7 | 0.685 6 | 0.654 5 | 0.606 1 | |
| VP | 0.703 4 | 0.690 2 | 0.660 5 | 0.611 1 |
从表 2中可以看出,从基于视觉感知的角度,本文算法无论是单项还是综合指标上,都优于小波变换方法、图像融合方法和图像修复方法,进一步说明了本文方法在去云效果上的优越性。
5 结束语由于采用结构风险最小化原理,支持向量机具有抗噪性能强、泛化性好等优点,在数据不确定、小样本下更具优势[15]。本文利用支持向量机进行云检测,并结合太阳角度信息检测云阴影;接着使用支持向量值轮廓波变换,通过二值图构造选择矩阵进行影像镶嵌的方法先完成厚云和阴影的初步去除。经过云层和云阴影初去除的图像中,还有小范围的云层及云阴影重叠区域,该区域不能通过影像镶嵌方法进行云层及云阴影去除,通过统计学补偿的方法可以进行修复。仿真试验表明,本文提出的基于支持向量机的云去除方法的视觉效果比小波变换方法、图像融合方法和图像修复方法要光滑且很清晰,并且定量评价标准也说明了本文方法在去云效果上的优越性。
| [1] | ZHU Xifang,WU Feng.An Improved Approach to Remove Cloud and Mist from Remote Sensing Images Based on Mallat Algorithm[C]//Proceedings of SPIE.Beijing:[s.n.],2007. |
| [2] | YANG J,ZHAO Z,MA J,et al.Image Fusion for Automatic Detection and Removal of Clouds and Their Shadows[C]//Proceedings of SPIE.San Jose:[s.n.],2006. |
| [3] | GABARDA S,CRISTOBAL G.Cloud Covering Denoising through Image Fusion[J].Image and Vision Computing.2007,25:523-530. |
| [4] | YANG Jun,ZHAO Zhongming,YANG Jian.A Shadow Removal Method for High Resolution Remote Sensing Image[J].Geomatics and Information Science of Wuhan University.2008,33(1):17-20.(杨俊,赵忠明,杨健.一种高分辨率遥感影像阴影去除方法[J].武汉大学学报:信息科学版,2008,33(1):17-20.) |
| [5] | MAALOUF A,CARRE P,AUGEREAU B,et al.A Bandelet-based Inpainting Technique for Clouds Removal from Remotely Sensed Images[J].IEEE Transactions on Geoscience and Remote Sensing,2009,47(7):2363-2371. |
| [6] | EL-ARABY E,EL-GHAZAWI T,LE MOIGNE J,et al,Reconfigurable Processing for Satellite On-board Automatic Cloud Cover Assessment[J].Journal of Real-time Image Processing,2009,4(2):245-259. |
| [7] | BENABDELKADER S,MELGANI F.Contextual Spatiospectral Postreconstruction of Cloud-contaminated Images[J].IEEE Geoscience and Remote Sensing Letters,2008,5(2):204-208. |
| [8] | ZHENG S,SHI W,LIU J,et al.Remote Sensing Image Fusion Using Multiscale Mapped LS-SVM[J].IEEE Transactions on Geoscience and Remote Sensing,2008,46(5):1313-1322. |
| [9] | HU Gensheng,LIANG Dong,KONG Jie.Remote Sensing Image Fusion Based on Support Vector Value Contourlet Transform[J].Acta Electronica Sinica,2010,38(6):1287-1292.(胡根生,梁栋,孔颉.基于支持向量值轮廓波变换的遥感影像融合[J].电子学报,2010,38(6):1287-1292.) |
| [10] | KONG Jie.Research of Cloud Removing Method for Remote Sensing Image Based on Support Vector Machine[D].Hefei:Anhui University,2010.(孔颉.基于支持向量机的遥感图像云层去除方法研究[D].合肥:安徽大学,2010.) |
| [11] | LI Xiaochun,WANG Yong,CHEN Jing.Detection and Removal of Clouds and Their Shadows from Multi-spectral Image[J].Journal of Astronautics,2004,25(5):555-559.(李小春,王勇,陈鲸.多光谱图像中云层及阴影的检测与消除[J].宇航学报,2004,25(5):555-559.) |
| [12] | TAN Kun,DU Peijun.Wavelet Support Vector Machines Based on Reproducing Kernel Hilbert Space for Hyperspectral Remote Sensing Image Classification[J].Acta Geodaetica et Cartographica Sinica,2011,40(2):142-147.(谭琨,杜培军.基于再生核Hilbert空间小波核函数支持向量机的高光谱遥感影像分类[J].测绘学报,2011,40(2):142-147.) |
| [13] | TSENG D C,TSENG H T,CHIEN C L.Automatic Cloud Removal from Multi-temporal SPOT Images[J].Applied Mathematics and Computation,2008,205:584-600. |
| [14] | ZHU Shulong,ZHU Baoshan,WANG Hongwei.Process and Application of Remote Sensing Image[M].Beijing:Science Press,2006:95-111.(朱述龙,朱宝山,王红卫.遥感图像处理与应用[M].北京:科学出版社,2006:95-111.) |
| [15] | WU Xin,ZHANG Huanlong,SHU Yunxing.Mosaic Image Quality Evaluation Method Based on Visual Perception[J].Computer Engineering,2008,34(18):220-222.(武新,张焕龙,舒云星.基于视觉感知的镶嵌图像质量评价方法[J].计算机工程,2008,34(18):220-222. |