1、利用 L 波段星载重复轨道干涉 SAR 提取DEM 及大气效应分析 薛东剑 郑洁 李成绕 李婉秋 苏璐璐 成都理工大学地球科学学院 中国林业科学研究院资源信息研究所 国土资源部地学空间信息技术重点实验室 摘 要： 以四川省茂县、安县地区的 ALOS PALSAR L 波段雷达影像为数据源, 在对地形复杂区的配准、解缠等算法研究的基础上, 结合 SAR 特有的成像几何结构, 对两幅 SAR 图像进行快速自动配准, 配准误差小于 0.2 个像元;在去除平地引起的干涉相位变化后, 运用 MCF 对缠绕相位进行解缠, 在此基础上提取研究区的数字高程模型 (DEM) , 分析发现精度受多种因素影响, 其
2、中波长较长的 L 波段数据比波长较短的时间相关性好, 此外为消除大气效应, 采用改进的相位累积法去除大气的影响, 在一定程度上提高 DEM 精度。其误差范围在10m 左右, 对快速提取地形信息具有一定的借鉴意义。关键词： PALSAR; InSAR; 配准; 相位解缠; 大气效应; 作者简介：薛东剑 (1977-) , 男, 博士.收稿日期：2016-11-23基金：四川省教育厅重点资助项目 (16ZA0100) Using repeat orbit L-band interferometric SAR to extract DEM and analyzing the atmospheric
3、effectsXUE Dongjian ZHENG Jie LI Chengrao LI Wanqiu SU Lulu College of Earth Sciences of CDUT; Abstract： This paper, taking ALOS PALSAR L-band radar images as the data source in the counties of Maoxian and Anxian, Sichuan Province, adopts a coherent coefficient method aiming at two SAR images automa
4、tic registration and MCF algorithm at phase unwrapping based on the research of registration, unwrapping algorithms, extracted study areas Digital Elevation Model (DEM) .The key technology of interferometric synthetic aperture radar is focused on the basis of the analysis of the current research sta
5、tus.At the time of interferometry, the SAR data with longer wavelength L-band has better temporal coherence than the shorter wavelength.Then the accuracy and characteristics of DEM data are analyzed, whose most error range is about 10 m, and the influence of atmospheric effects in interferometriv ap
6、erture radar topographic map is discussed, which can provide the basic data for the research work of subsequent surface deformation.Keyword： PALSAR; InSAR; image registration; phase unwrapping; atmospheric effect; Received： 2016-11-23以往运用遥感技术提取地形信息, 多采用光学图像, 通过立体相对进行提取;图像畸变及天气条件严重限制了其应用。到 20 世纪 60 年
7、代末, 合成孔径雷达干涉测量技术 (InSAR) 应运而生, 其能够全天候、全天时工作, 精确度高, 具有光学图像无法比拟的优势, 有着广阔的应用前景。Graham (1974) 等1利用机载 SAR 数据获取满足 125 万地形图要求的高程数据, 开创 InSAR 技术在对地观测中获取三维信息的先河。随着星载合成孔径雷达的不断升空, 产生大量丰富的 SAR 数据, 雷达干涉被广泛地应用于提取地形信息2-5。研究表明, 波长较长的 L 波段 SAR 数据比波长较短的 X、C 波段数据, 在干涉测量时, 时间相关性好, 更适合于地形起伏大、植被覆盖密、形变梯度大的区域3, 尤其有利于西南地区高植
8、被覆盖区地形提取。并且 ALOS PALSAR FBS 数据的干涉临界基线距较大, 能够使用更多的影像参与干涉处理。由于 SAR 特有的成像几何结构, 其处理中的很多关键技术还一直处于不断研究与发展中, 尤其针对不同的研究区, 其影响提取精度的处理是目前研究的焦点与热点 (如解缠、大气效应等) 。采用茂县地区的 L 波段 ALOS PALSAR 数据, 在对配准、解缠算法等研究的基础上, 提取了数字高程模型 (DEM) , 并重点分析了大气效应对提取精度的影响。采用该技术可快速提取地形信息, 这对于地形变化监测及由地形引起的灾害监测预报、环境治理等问题具有重要的实用价值。1 研究区概况及数据预
9、处理研究区位于四川省东北部的茂县、安县地区, 地处龙门山中南段与四川盆地相结合区域, 地理坐标:东经 10335471042240, 北纬 312348320013, 是“5.12”地震中受灾较为严重地带, 对其震后地形的提取具有一定的应用价值。研究区地形主要为平坝、丘陵 (台地) 、中低山三种类型, 既以大光包斜冲断层和北川冲断层为界, 西北部属四川西部地槽区的前龙门山褶断带, 系龙门山脉, 地势较高, 最高峰位于九顶山主峰高 4 989m, 为区内最高点, 最低点位于安县秀水镇, 海拔为 548m。研究中采用 2010-01-24 与 2010-03-11 1.1 级 L 波段 (波长 :
10、0.236m) ALOS PALSAR 数据, 每景数据覆盖范围约为 2 565.2km。主辅单视复数 (SLC) 影像距离向像元大小 4.68 m, 方位向 3.18m, 9 344 列, 18 432 行 (见表 1) , 其中以 2010-01-24 数据为基准 (升轨) , 进行干涉处理 (见图 1, 为研究区 ALOS PALSAR HH 数据, 为显示方便, 在此进行了多视处理及地理编码) ;参考 DEM 数据选择美国 NASA 和 NIMA 联合测量的 SRTM 数据。图 1 研究区 ALOS PALSAR 数据示意 下载原图表 1 ALOS PALSAR 主辅影像参数 下载原表
11、 2 InSAR 关键处理技术利用星载 SAR 重复轨道方式进行干涉测量, 主要是通过重复两次观测, 获得同一区域的单视复数影像对, 根据两幅天线和观测目标之间的几何关系 (见图 2) , 结合观测平台的轨道参数, 在对相位差解算的基础上, 提取地面高程。由图中几何关系可推算高程:图 2 雷达干涉测量原理示意图 下载原图图 2 中, A 1, A2为天线在两个成像时刻的位置;H 是 A1的高度;B 是基线, 为基线与水平方向的夹角; 为天线在 A1处的视角, 为波长;, + 为 A1, A2到观测点 Z 的斜距;h 为目标点高程。利用重复轨道获取的相位差仅是相位差值中的小数部分, 需联合相位解
12、缠求出 t。由于 SAR 特有的成像方式, 使其在数据获取上具有一定要求, 且处理上也有一定的难度, 尤其是对两幅图像的配准及相位解缠的好坏直接关系到提取地形的精度。2.1 SAR 图像自动配准对不同天线获得数据进行精确配准, 估测两幅 SLC 影像之间距离和方位向的偏差, 是利用 SAR 进行干涉测量提取 DEM 的关键步骤。为了准确计算对应的地面目标的干涉相位或形变相位, 需要实现雷达图像之间亚像元级精度的配准12。在光学图像上, 自动配准的基本思想是相同的地物具有相似的光谱特征, 通过对两个图像做相对移动, 找出其相似性量度值最大或差别最小的位置作为图像配准的位置。SAR 图像的配准也是
13、类似的原理, 目前国内外大部分配准方法基本上是在确定控制的基础上, 以控制点坐标拟合多项式:其中, x r, yr为参考影像中控制点坐标;X t, Yt为待配准数据相应控制点坐标。如何精准确定控制点, 是配准处理的前提, 也是关键, 文中针对 ALOS PALSAR成像特征, 综合分析相位差平均强度梯度法、相干系数法、绝对差值法、最大频谱法的基础上, 采用相干系数法进行了控制点选取, 此方法充分利用了振幅与相位信息, 在一定的搜索范围内, 逐一计算每一点的相关系数, 以最大值为配准点的位置。为了避免产生模糊问题, 取得精确的偏差估计, 在多视估计后, 进行一次单视处理, 并采用改进的最小二乘法
14、生成配准多项式 (表 2) , 从而使配准误差小于 0.2 个像元, 提高干涉相关性。表 2 干涉相对配准拟合参数 下载原表 2.2 干涉图生成及相位解缠对干涉像对进行配准后, SLC1 与 SLC2 影像的数据可表示为通过共轭相乘, 可以计算出同名点上的相位差:其中, 表示复数的共轭, int的相位是每一同名点上相应于两个观测点的相位差 () , 所得相位的主值数据就是被缠绕相位数据, 以干涉图形式表示。运用 InSAR 技术提取 DEM 的精度是否成功, 关键是干涉像对的相位相关性的高低。经处理, 在本次实验中, 相干系数最大值是 0.989, 最小值是 0.001 4, 平均值是 0.6
15、91 6, 0.4 以上相关性占 78%, 该试验区干涉像对的整体相干性较好 (图 3) 。此时生成的干涉图中是由多个成分组成, 需去平处理, 把干涉图中的地球弯曲距离和方位向的相位趋势去除, 去除了平地相位之后, 干涉条纹由密集变为稀疏, 避免了将平地相位信息作为误差引入高程和形变信息。从图 4 可以看出, 去平后干涉图, 在地形比较平坦的地区 (图 4 右上角地区) 干涉条纹已经去除, 这与实际地形比较相符。图 3 相干图与强度合成 下载原图图 4 去平后干涉 下载原图相位解缠在干涉处理中具有重要的地位, 而方法的选择会很大程度影响高程的精度。常用的算法有区域增长法, 最小费用流量法等,
16、其中, 区域增长算法, 是根据相位的质量来规划积分路径, 先选择高质量的区域进行积分, 然后再向低质量的区域, 从而得到解缠图像。这种方法在相位图质量较高时, 解缠结果好, 速度快, 但无法识别残差点, 会引入相位误差影响精度。最小费用流法是一种基于网络规划的相位解缠方法, 其基本思想是最小化解缠后相位的导数与缠绕相位的导数之间的差异, 该方法可兼顾精度与处理效率, 尤其对范围较大地区的解缠具有一定的优势。在此以最小费用流量 (MCF) 算法设置的阀值作为标准 (在此相干阈值设为 0.3, 图 5) , 将所有相干性满足条件的相位添加到一个集合, 在集合中建立三角网, 构造对偶图, 用最小费用
17、流法连接各残差点, 求出最小费用流集合, 再对相位图进行积分, 根据相位质量从高到低进行解缠。图 5 相位解缠结果 下载原图3 地形信息提取及大气效应分析经过相位解缠图像 (图 5) 的相位主值已经转化成了真实的相位, 采用式 3 计算出每个像元的高程, 在此基础上将获取的高程信息由雷达坐标转化为地理坐标 (图 6) 。为了分析提取出的 DEM 精度, 采用 SRTM DEM 数据做参考, 进行比较分析。SRTM 数据在平原区精度优于 1250 000DME, 中误差为 4.921。运用 InSAR 技术提取的 DEM 与 SRTM DEM 数据比较发现, 整体趋势比较吻合 (图 6) , 平
18、坦地区误差范围较小, 而坡度较大的地区, 误差相对较大, 存在部分异常值, 但差值大的像元比较少。 (1) 研究区大部位于山区, 有植被覆盖的山区受时间失相关的影响, 对提取的精度造成一定的影响, 而城市等平坦地区受时间失相关的影响较小, 其相关性较好; (2) 雷达干涉测量中, 除了时间去相干外, 大气效应也是重要的影响因素, 不同时间获取的图像受大气影响产生不同的相位延迟, 干涉相对中两个 SAR 图像回波信号的相位可表示为图 6 InSAR 提取的 DEM 与 SRTM DEM 叠合分析 下载原图其中, 1, 2为斜距, 1, 2为雷达电磁波穿过大气层 (电离层与对流层) , 受电子浓度
19、、气压、温度、水汽等影响造成的大气延迟, 为波长。 为两次观测大气引起的相位差:大气造成的延迟主要发生在电离层与对流层。以往研究表明, 大气影响主要是由于水汽的变化引起的17-18, 表征大气水汽的物理量主要为可降水汽含量PWV (Precipitable Water Vapor) , 可以转换为天顶湿度延迟 ZWD (Zenith Wet Delay) :其中, K 为折射常数, w为液态水密度, T M为对流层的加权平均温度, R v为大气常数。针对重复轨道干涉测量中高程与相位之间的参量关系, 可推导出: h为高程误差, 为天顶湿度延迟误差, 为入射角, B 为垂直基线。综合目前研究现状1
20、5-20, InSAR 测量中常用的减少大气影响的方法, 除了数据本身因素外, 主要有相位累积法, PS 技术, 基于 GPS 数据的改正方法以及利用外部水汽产品 (如 MERIS、FY、MODIS 等) 校正法等, 其中 PS 技术、GPS 校正及外部水汽校正法, 虽然普适性较高, 但主要适用于地表形变监测中21。结合 ALOS PALSAR 数据特点, 且考虑到外部数据获取的同步性及相位解缠的难易等因素, 文中借鉴相位累积法, 采用 2010 年 ALOS PALSAR 数据生成三对独立干涉图进行相位梯度累加分析, 相位梯度可表示为R 和 I 是去除了平地效应后干涉图的实部和虚部, 垂直基
21、线为 B 的 N 个干涉图的平均相位梯度 表示为为减少噪声的影响, 对干涉图的相位梯度减去平均值和垂直基线的乘积, 得到相位梯度变化则大气引起的相位延迟异常 为为相位梯度异常, i为空间相关系数。通过该相位梯度累积平均一定程度上降低了大气的影响。大气改正后的高程差异值明显低于未经过大气改成的高程差异值, 且高程差异分布相对比较集中, 高程差异的标准差由15.65 m 提高到了 6.02m, 尤其是平坦地区误差范围在10m 以内的置信度由68.78%提高到了 93.21%, DEM 精度得到提高。4 结束语以 ALOS PALSAR 数据为例, 在对 InSAR 处理流程及算法详细分析的基础上,
22、 研究了地形复杂区域 SAR 图像配准、解缠等核心算法;以最小费用流 (MCF) 算法为基础, 将相干性满足阈值条件的相位添加到一个集合, 构造对偶图, 用最小费用流法连接各残差点, 对相位进行积分, 根据相位质量从高到低进行解缠, 并提取了 DEM。研究表明利用雷达干涉测量提取地形信息具有自动化程度高、效率较高等特点, 但失相关是限制 INSAR 技术应用的一个较为严重的问题, 尤其植被覆盖度大、地形复杂地区, 时间去相关明显;除此外大气效应也是重要的影响因素, 在地形提取中, 运用相位梯度累加在一定程度上可以消除大气的影响, 提高 DEM 的精度;ALOS 数据由于波长较长, 对地表穿透深
23、度大, 干涉性相对较好, 有利于西南地区高植被覆盖区地形提取。参考文献1GRAHAM L C.Synthetic Interferometer Radar for Topographic MappingA.Proceeding of the IEEE, 1974, 62 (6) C.s.l.:s.n., 1974.763-770. 2GOLDSTEIN R M, ZERBER H A, WERNER C L.Satellite Radar Interferometry:Twodimensional PhaseUnwrappingJ.Radio Science, 1988, 23 (4) :71
24、3-720. 3RAUCOULES D, COLESANTI C, CARNEC C.Use of SAR interfer ometr y fo r detecting and assessing g ro und subsidenceJ.C.R.Geo science, 2007 (339) :289-302. 4ZEBKER H A, ROSEN P A, Hensley S.Atmospheric Effects in Interferometric Synthetic Aperture Radar Surface Deformation and Topographic MapsJ.J
25、ournal of Geophysical Research, 1997, 103 (B4) :7547-7563. 5HSIEH Chia-Sheng, SHIH Tian-Yuan, CHING Jyr, et al.Using differential SAR interferometry to map land subsidence:a case study in the Pingtung Plain of SW TaiwanJ, Nat Hazards (2011) 58:1311-1332. 6KOCH A, HEIPKE C.Quality assessment of digit
26、al surface models derived from the Shuttle Radar Topography Mission (SRTM) A.the IEEE 2001International Geoscience and Remote Sensing SymposiumC.Sydney, Australia:2001. 7DEO R, MANICKAM S, RAO Y S, et al.Evaluation of interferometric SAR DEMS generated using TanDEMX data.2013IEEE International Geosc
27、ience and Remote Sensing Symposium (IGARSS) C.Melbourne, VIC:IEEE, 2013:2079-2082. 8陈尔学, 李增元, 车学俭.SAR 干涉测量数字高程模型提取与高程误差校正J.高技术通讯, 2000 (7) :57-63. 9刘国祥, 丁晓利, 李志林, 等.InSAR 建立 DEM 的试验研究J.测绘学报, 2001, 30 (4) :337-343. 10吴宏安, 张永红.ALOS/PALSAR InSAR 数学模型研究J.遥感信息, 2011 (2) :106-109. 11沈强, 乔学军, 金银龙, 等.ALOS P
28、ALSAR 雷达影像 InSAR 数据处理中的基线和地形误差分析J.大地测量与地球动力学, 2012, 32 (2) :1-6. 12刘广, 郭华东, 范景辉.基于外部 DEM 的 InSAR 图像配准方法研究J.遥感技术与应用, 2008, 23 (1) :72-76. 13LEE J S.Digital Image Smoothing and the Sigma Filter.ComputerVisionJ.Graphies and Image Proeessing, 1983, 24:255-269. 14WANG Qi-jie.Generalized functional model
29、of maximum and minimum detectable deformation gradient for PALSAR interferometryJ.Trans.Nonferrous Met.Soc.China 24 (2014) 824-832. 15HANSSEN R F.Radar Interferometry:Data Interpretation and Error AnalysisM.Dordreeht:Kluwer Academic Press, 2001. 16徐佳, 关泽群, 何秀凤.重复轨道雷达干涉测量中的大气影响及其研究展望J.国土资源遥感, 2007, 7
30、2 (2) :1-7. 17ZEBKER H A, ROSEN P A, HENSLEY S.Atmospheric effects in interferometric synthetic aperture radar surface deformation and topographic mapsJ.Journal of Geophysical Research, 1997, 102 (4) :7547-7563. 18孙广通, 张永红, 吴宏安.合成孔径雷达干涉测量大气改正研究综述J.遥感信息, 2011 (4) :111-116. 19SONG Xiaogang, LI Deren,
31、SHAN Xinjia, et al.Correction of atmospheric effect in ASAR interferogram using GPS and MODIS dataJ.Chinese Journal of Geophysics, 2009, 52 (6) :1457-1464. 20PUYSSEGUR B, MICHEL R, AVOUAC J P, et al.Tropospheric Phase Delay in Interferometric Synthetic Aperture Radar Estimated from Meterological Mod
32、el and Multispectral ImageryJ.Joural of Geophysical Research, 2007, 112:120-135. 21崔喜爱, 曾琪明, 童庆禧.重轨星载 InSAR 测量中的大气校正方法综述J.遥感技术与应用, 2014, 29 (1) :9-17. 22GOLDSTEIN R.Atmospheric limitations to repeattrack radar interferometry, G eophys R.es.L ett., 1995, 22 (18) :2517-2520. 23SANDWELL D T, PRICE E J.
33、Phase Gradient Approach to Stacking InterferogramsJ.Journal of Geophysical Research, 1998, 103 (B12) :30183-30204. 24FERRETTI A, PRATI C, ROCCA F.Multibaseline InSAR DEM Reconstruction:the Wavelet Appro achJ.IEEE Transactions on Geoscience and Remote Sensing, 1999, 37 (2) :705-715. 25Ramon Brcic, Alessandro Parizzi, Michael Eineder, et al.Estimation and compensation of onospheric delay for SAR interferometryC.IEEE International Geoscience&Remote Sensing Symposium, 2010, 38 (1) :2908-2911.
