王睿,韦春桃,马云栋,胡涛
(1.桂林理工大学测绘地理信息学院,广西桂林541006;2.广西百色开发投资集团有限公司,广西百色533000;3.重庆交通大学土木建筑学院,重庆400074)
光学影像的波长设置一般较短,地物的光谱信息易受天气影响,常被云层遮挡,这在很大程度上缩小了光学遥感影像的选择范围。因此,探讨如何修复带云影像在遥感应用方面具有重要意义。
带云影像的修复根据云层薄厚来选择处理方法:薄云区影像由于存在云层和下垫面的混合像元,常采用滤波法[1]去除处于低频的云层信息来修复影像;厚云区影像完全遮挡下垫面信息,需参考多时相影像进行地物的还原,修复的主要步骤包括云区、阴影的检测以及去除。常用的阴影检测方法是通过云增强模型来突出云区光谱信息,再设定阈值检测云:掩膜方程法(Fmask)[2]需要云层高度数据且计算量大;色彩空间变换法[3]依赖热红外波段,通用性较弱;多时相云检测法(MTCD)[4]需要较为完整的影像数据建立厚云增强模型,且易受地物变化影响。在去除方面主要是依靠多时相遥感影像进行替换,多数研究采用线性模型解决影像间的光谱差异的问题[5-6]。考虑到大气条件、季节变化、地类反射差异等因素影响,影像间的光谱差异会在局部区域偏离线性变化规律,导致整体影像存在非线性关系。
针对多时相遥感影像之间光谱信息存在非线性问题,本文采用BP神经网络算法拟合两幅影像之间光谱变换的映射关系,将无云参考影像与带云目标影像的光谱信息进行匹配,减少影像间的光谱差异,再把光谱纠正后的参考影像镶嵌入云和阴影的检测区域,修复含有厚云的遥感影像。
本文去除Landsat影像上的厚云及其阴影方法可分为3个步骤:①利用HOT算法检测云区、阴影成像几何关系检测阴影,得到云区与阴影的检测图;②设定样本和期望值训练BP神经网络,将无云参考影像输入网络后得到与带云目标影像光谱信息相近的图像;③将光谱匹配后的图像对应云和阴影检测的像元位置镶嵌入带云目标图像,最终去除云和阴影。
Landsat数据中蓝波段与红波段相关性高,在晴空范围内两者的二维散点图呈现出一条带状的“晴空线”;但有云层的影像成像条件不佳,各地物在蓝波段处的反射率比红波段处的反射率升高的更快,使晴空线上的地物点向蓝波段轴靠近。厚云在各个可见光波段处反射率均较高、反射变化程度小,云像元集中在另一条带状直线上,与晴空线形成一个夹角,这里定义它为“云层线”,此时散点图形成“箭头”状(图1a)。根据云层线偏离晴空线的特点,从其物理模型分析偏移量(图1b)可以将厚云与其他地物区别开,这就是haze optimized transform[7],简称HOT算法,计算公式为[8]
式中:Bblue、Bred分别是地物像元在蓝色和红色波段的灰度值;θ为晴空线倾角;AK为晴空线的截距。
带云影像在居民地处常有与云层反射率相似的高亮体,高亮地物偏离晴空线的量值小于云层偏移量,使“高亮体线”位于晴空线与云层线之间。导致HOT算法在计算云层偏移量时涵盖了高亮体偏移量,使检测结果中二者混淆。本文将晴空线倾角θ换为高亮体线倾角,以高亮体线为基准计算HOT(即云层线到高亮体线的距离)。这样就避开了高亮体的干扰,得到较为准确的云区检测图。
图1 云检测示意图Fig.1 Diagram of cloud detection
厚云与云影常成对出现,云影相对于云的位置偏向基本保持一致,并且各个云影伸缩变化的幅度范围类似。检测云影时在目标影像上选择n个分布均匀的碎云Mi(i=1,2,…,n),获取其多边形质心坐标:碎云Mi(x1i,y1i)与其对应的阴影质心坐标Ni(x2i,y2i)(图2a),计算两个质心之间的距离di和连线MiNi的方位角α,统计n对质心间距的平均值和标准差σ,则对于整幅影像来说阴影距离厚云的偏移量为d=+3σ[9]。在云区检测斑块图的基础上加上偏移量进行偏移,并采用数学形态学膨胀算法解决阴影形变的问题(图2b),最终得到阴影检测图。
图2 阴影检测示意图Fig.2 Diagram of shadow detection
人工神经网络是计算机对人脑或自然神经网络若干个基本特征的抽象和模拟,把大量简单的处理单元连接成网络结构,解决高度复杂的非线性动力学问题[10-11]。其中在遥感应用领域得到广泛使用的神经网络是误差反向传播学习网络(Back Propagation network,简称BP网络)。
1.3.1 光谱匹配原理分析及操作BP网络算法首先需要人工获取样本值P和期望值T,通过网络激励函数的计算将样本值P向期望值T转换。以输出误差为依据调整网络权值和阈值,使网络不断逼近样本值与期望值之间的映射关系。
为减少计算量,本文以影像灰度值代表其光谱特性并参与计算。BP算法在光谱匹配实验中的实质是把求解两幅遥感影像的光谱特性转换问题转化为求解样本灰度值与期望灰度值之间的转换关系问题。当获取晴空区无云参考影像和带云目标影像的光谱特征点时,该点所表示的灰度值分别作为样本值P和期望值T输入BP网络。样本值P经过正向传播计算验证输出值与期望值T的均方误差是否超限,再由误差反向传播优化网络参数,往返计算,直到误差达标或是循环次数达到设限为止。此时,训练好的网络已保存各个计算层最优的权值和阈值,最终获取辐射特征之间的转换关系。光谱匹配的具体操作如下:
①将参考影像和目标影像进行配准,保证选取特征点时两幅影像在同一坐标点下地类相同;
②在参考影像上进行ISODATA聚类,以判别的地物类型为依据,选取晴空区各地类纹理变化明显处作为特征点,则样本值P由待修正光谱信息的参考影像特征点灰度值组成,期望值T由样本特征点对应地类在目标影像上的灰度值组成;
③训练BP网络后获取光谱转换关系,将整幅无云参考影像输入网络,根据优化的网络参数计算可得到与目标影像光谱相近的图像。
1.3.2 网络训练结构为减少计算量及缩短计算时间,本文采用单个隐含层,可基本满足实验要求。实验中使用参考影像和目标影像各自对应的3个波段分别作为输入、输出数据,因此输入和输出节点数均设为3个。隐含层节点数的选择目前尚无通用的理论指导,经反复实验及调试认为10个节点对实验效果较为理想。因此利用BP神经网络进行光谱匹配的网络结构设计为3-10-3。
无云参考影像经过BP神经网络光谱匹配后,其光谱特征已与带云目标影像相近。对应云区和阴影检测图相应的像元位置替换目标影像的像元去除云区和阴影。
为验证参考影像与目标影像的成像时差大小对实验结果是否有影响,本文采用两组Landsat 8数据的波段7、5、3影像(大小均为1 024×1 024像元,分辨率30 m)计算,以说明算法的有效性。
实验1。使用桂林部分地区2013-09-01和2013-09-17两幅影像,图3a、3b分别给出3个波段的假彩合成图。该区域地表覆盖物主要为植被、水体、裸地、居民地和田地,由于时相相近,两幅影像的地物类型未发生明显变化。目标影像中厚云与阴影的下垫面均包含这几个地类,以此检验实验结果对不同纹理区域修复的适用性。图3c中白色区域为厚云检测位置,灰色区域为阴影检测位置,图3d为目标影像的最终修复结果。可见,BP神经网络算法在裸地、居民地和田地处的细节表达明显,修复图像较为清晰。
实验2。使用南宁部分地区2013-09-15和2013-11-18两个Landsat 8数据,图4a、4b为3个波段的假彩合成图。由于两幅影像的成像时间存在季节变化,植被和耕地地区发生较大改变。经过BP神经网络的光谱信息校正,能比较准确地预测出被厚云和阴影遮挡地区(图4c为检测结果)原本的地物像元灰度值,对目标影像进行镶嵌后图像得到了较好的改善,如图4d所示。
图3 地物未发生变化去云结果Fig.3 Unchanged object in cloud detection
图4 地物发生变化的去云实验结果Fig.4 Changed object in cloud detection
目标影像上的厚云和阴影区域是通过像元替换镶嵌进行修复的,本文在修复结果图上抽取云区的镶嵌部分像元与晴空区像元,计算其灰度均值、信息熵和平均梯度。数值越接近目标影像值,则说明镶嵌区像元的色彩表现程度、细节、表现力和清晰度越好。
从图像整体的镶嵌质量方面看,根据空间频率活动性(SFA)[12],分别按照修复后图像的行、列计算灰度值I(i,j)之间的差异,以反映镶嵌边界的拼缝是否光滑。SFA值越小,表明在接缝处的像元灰度值过度得越平缓、拼接质量越高。其公式定义为
表1是对比本文BP算法与文献[6]的线性模型算法在波段7处去云方面的效果,可看出BP算法在均值、信息熵和平均梯度方面更接近目标影像,而SFA值则小于线性模型,说明神经网络的镶嵌边界更为平滑,在去除厚云及其阴影方面具有一定的可行性和有效性。
表1 波段7厚云及其阴影镶嵌结果评价Table 1 Comparison of cloud and shadow removing methods in Band 7
利用BP神经网络进行图像的光谱匹配,对参考影像颜色校正后能保持原有分辨率,光谱特征得到统一、纹理清晰。将高亮体线倾角替换晴空线倾角计算HOT,能够在一定程度上避开高亮体对云检测的干扰。利用云区检测图进行偏移和膨胀计算检测云影,无需再确定阈值进行分割,算法简单有效。
由于大幅影像的云影偏移差异较大,需要分区域处理且人工干预较多,仍需要进一步的研究。
[1]李刚,杨武年,翁韬.一种基于同态滤波的遥感图像薄云去除算法[J].测绘科学,2007,32(3):47-48.
[2]Zhu Z,Woodcock C E.Object-based cloud and cloud shadow detection in Landsat imagery[J].Remote Sensing of Environment,2012,118:83-94.
[3]李微,李德仁.基于HSV色彩空间的MODIS云检测算法研究[J].中国图象图形学报,2011,16(9):1696-1701.
[4]Hagolle O,Huc M,Pascual D V,et al.A multi-temporal methodforclouddetection,appliedtoFORMOSAT-2,VENμS,LANDSAT and SENTINEL-2 images[J].Remote Sensing of Environment,2010,114:1747-1755.
[5]周伟,关键,姜涛,等.多光谱遥感影像中云影区域的检测与修复[J].遥感学报,2012,16(1):132-142.
[6]李炳燮,马张宝,齐清文,等.Landsat TM遥感影像中厚云和阴影去除[J].遥感学报,2010,14(3):534-545.
[7]Zhang Y,Guindon B,Cihlar J.An image transform to characterize and compensate for spatial variations in thin cloud contamination of Landsat images[J].Remote Sensing of Envi-ronment,2002,82(2-3):173-187.
[8]李存军,刘良云,王纪华,等.基于Landsat影像自身特征的薄云自动探测和去除[J].浙江大学学报:工学版,2006,40(1):10-13.
[9]秦雁,邓孺孺,何颖清,等.基于光谱及几何信息的TM图像厚云去除算法[J].国土资源遥感,2012(4):55-61.
[10]Bao Y H,Ren J B.Wetland landscape classification based on BP neutral network in DaLinor Lake Area[J].Procedia Environmental Sciences,2011,10:2360-2366.
[11]Zhang L P,Wu K,Zhong Y F,et al.A new sub-pixel mapping algorithm based on BP neural network with an observation model[J].Neurocomputing,2008,71(10-12):2046-2054.
[12]武新,张焕龙,舒云星.基于视觉感知的镶嵌图像质量评估[J].计算机工程,2008,34(18):220-222.