基于NEST的差分干涉测量及其应用

杜顺季，邢 诚,2

（1.武汉大学 测绘学院，湖北 武汉 430079；2.精密工程与工业测量国家地理信息局重点实验室，湖北 武汉 430079）

合成孔径雷达作为主动式微波遥感技术，因不受时间和天气条件的限制，成为对地观测系统重要组成部分。1993年Massonnet等人在《Nature》上发表了利用ERS-1获得与GPS一致的Landers地震形变场的文章[1]，由此DInSAR技术得到了广泛研究和应用。

目前能够系统处理SAR数据软件的主要是商用Gamma、ENVI的SARscape模 块 和doris、ROI_PAC等开源软件，商用软件价格不菲而开源软件一般基于Linux平台的命令模式，不方便用户操作和交互。本文利用以欧空局为主研发的跨平台开源软件NEST来进行SAR图像处理的研究，以伊朗Bam地震为例，利用ASAR影像基于两轨法DInSAR技术获取视线向的地表形变量，介绍了NEST软件处理SAR数据的过程，得到了较好的结果。

1 二轨法DInSAR基本原理

根据去地形影响所采用的数据和方法，DInSAR可分为二轨法、三轨法和四轨法[2]。二轨法是利用外部高精度地面数字高程模型(DEM)和轨道数据来模拟地形相位，从原始干涉图中减去该模拟的地形相位就得到新的差分干涉相位图，去除了地形影响而仅包含了形变相位。

利用二轨法进行形变相位的提取处理比较简单，其优点是无需对干涉图相位解缠，避免了解缠困难[3]。在没有外部DEM的情况下，可选择三轨法和四轨法。但是利用三轨法和四轨法的地形像对获取的DEM会进一步引入大气效应的影响，从而引起形变测量的误差。二轨法DEM高程误差会直接影响形变相位，根据SAR测量几何图形和误差传播率可知：

式中，σΔr为形变测量误差；σh为DEM高程误差；B⊥为垂直基线分量；R为地面点至天线距离。由式(1)可知，当垂直基线为0时，DEM高程误差对形变测量结果无影响，这在地基SAR中是重要特性[4]，可以直接利用2幅图像获得地表形变，不需要考虑测量区域的地形相位，但是在星载SAR中基线不可能为零。目前美国航空局公开了全球的SRTM数字高程数据，高程精度可达16 m[5]，在平原地区绝对高程精度优于5 m，基本上满足了差分干涉测量的要求，为二轨法DInSAR的广泛应用提供了数据[6]。

2 NEST数据处理实践

2.1 NEST软件主要功能及特点

最新发布的NEST版本为4C-1.1，能够对SAR数据进行读取、后处理、分析和可视化，支持ERS1/2、Envisat的ESA卫星数据和诸如JERS、ALOS、TerreSAR-X、Radarsat1/2和Cosmo-Skymed等其他数据。主要功能包括有效的数据查看和管理、绝对定标、轨道与滤波处理、自动配准、干涉、地形改正、统计分析、海洋工具和全面集成doris。其特点为：良好的操作界面和灵活的参数设定、数据可视化和管理、自动下载轨道和DEM数据、批处理功能、数据统计和分析工具、多种数据格式输出等。该软件集成了doris软件算法，既保证了良好的图形化用户界面，又保留了高效的批处理框架，实现了SAR处理的各方面功能并提供众多分析工具，为SAR研究带来极大的方便[7]。

2.2 像对数据选取与处理

欧空局在Bam地震后公开了震区的Envisat ASAR数据，为VV极化，地震时间为2003-12-26，这次使用的2景降轨的地震前后图像，其相关参数见表1。

出该像对的垂直基线只有0.6 m，根据式(1)，用该数据计算形变场，地形高程引起的误差影响值很小，能够达到很高的形变监测精度。DInSAR数据处理的主要步骤包括：影像配准、干涉图生成、DEM模拟相位、干涉差分、干涉图相位滤波、相位解缠、形变量计算和地理编码。NEST软件实现了以上步骤并在每个操作时给出了详细的参数控制选项，能根据要求进行灵活的影像处理，相关的滤波和去噪处理可以选择性进行，其大致流程如图1所示。

图1 NEST处理流程图

数据输入阶段，NEST软件可通过Apply orbit file操作更新元数据中的轨道状态矢量。对于ASAR数据提供doris或delft精密轨道文件的自动下载，这些文件提供每60 s的卫星位置和速度。预处理中辐射校准是为了得到反映地物后向散射系数的图像，是关联像素值和散射系数的重要工作[8]，对于SAR数据的定量分析有重要意义。

在InSAR处理中，图像对的精确配准是至关重要的，保证2幅复图像的点对应地面同一个点，一般达到1/8像元的亚像元级配准精度才符合处理要求。NEST通过控制点选取操作进行精配准设置，可优化插值函数、配准窗口大小、多项式次数等参数，提供了最邻近、线性、4点/6点立方卷积、sinc插值等众多插值函数。本文采用精配准，二次多项式和6点立方卷积插值得到的配准残差均值Mean=0.02，标准差Std=0.02，达到了很高的配准精度，体现了该软件算法的准确性。

配准之后可对图像进行斑点噪声滤波，以便得到更好的干涉图。NEST提供2种滤波方式：针对单幅图像的经典滤波(包括中值滤波、均值滤波、Frost滤波、Lee滤波和Gamma滤波等)和多时相滤波，分别采用5×5的Gamma滤波和多时相滤波进行处理得到的图像如图2、图3所示。可以看出，滤波总体上有平滑的效果，特别是单幅图像的空间滤波使得图像更平滑但降低了分辨率，多时相滤波保留了一定的边缘信息，图像分辨率没有下降很多。多时相滤波引入时间维度，能够综合多幅图像中的信息，将相同位置不同时相的像元以某种方式进行组合，达到更好的斑点噪声抑制效果[9]。

图2 单幅图像滤波

图3 多时相滤波

经过配准的SAR图像共轭相乘得到干涉图，再消除平地相位，该相位的存在使得条纹密集，影响后续处理和分析，必须消除。NEST软件在同一个操作中实现，基于轨道参数和成像区域信息，通过分布在整幅影像的控制点进行最小二乘拟合计算平地相位，在干涉图中减去该部分相位。同时SAR成像的斜视几何效应，产生透视收缩、阴影、顶底倒置等图像畸变，二轨差分中NEST自动下载相应地区的DEM数据，依据轨道和影像数据进行地形改正，同时可完成地理编码，支持ENVI、HDF5、GeoTIFF等格式输出，特别是可输出为Google Earth支持的KMZ格式，方便用户直观分析结果。以上所有操作也可以在Graph Builder的命令行可视化工具中完成批处理，该窗口分为了2部分，上半部分是所执行命令过程图，下半部分是每条命令的详细参数设置，也是该软件强大的特色工具之一。

经过这一系列处理可得到二轨差分结果见图4，每条干涉条纹表示28 mm形变量，南面10个条纹表示卫星视线向形变约28 cm，北面7个干涉条纹，形变约19.6 cm，最大的相对位移变化大约为48 cm，与德国地学研究中心(GFZ)等机构公布的结果一致[10]，深入分析见文献[11]，输出到Google Earth叠加结果如图5所示。

图4 Bam形变干涉相位图

图5 地理编码输出到Google Earth结果图

3 结 语

本文介绍了二轨法差分干涉的主要步骤，并对Bam震区的SAR数据进行了处理。结果表明，DInSAR技术能有效获取地表变形，且NEST软件可以有效地进行SAR数据处理，良好的可视化操作界面使其使用更加方便，自动下载轨道和DEM数据并可对每个步骤详细设置，支持批处理，更加便捷和灵活，让SAR数据处理更加简单，有利于相关问题的学习和研究。

[1]Massonnet D, Rossi M, Carmona C, et al. The Displacement Field of the Landers Earthquake Mapped by Radar Interferometry[J]. Nature,1993, 364(6433): 138-142

[2]廖明生，林珲. 雷达干涉测量: 原理与信号处理基础[M]. 北京: 测绘出版社, 2003

[3]张红，王超，吴涛，等. 基于相干目标的DInSAR方法研究[M]. 北京: 科学出版社, 2009

[4]Casagli N, Catani F, Del Ventisette C, et al. Monitoring,Prediction, and Early Warning Using Ground-based Radar Interferometry[J]. Landslides, 2010, 7(3): 291-301

[5]陈俊勇. 对SRTM3和GTOPO30地形数据质量的评估[J].武汉大学学报:信息科学版,2005, 30(11): 941-944

[6]何敏，陆晓燕，何秀凤. 利用DInSAR二轨法监测徐州大屯中心区地表形变[J]. 地理空间信息,2011,9(5): 3-5

[7]ESA. NEST User Manual .[EB/OL].http://nest.array.ca/web/nest/documentation,2012-02-21

[8]袁礼海，葛家龙，江凯,等. SAR辐射定标精度设计与分析[J]. 雷达科学与技术,2009, 7(1): 35-39

[9]龚丽霞，吴樊，曾琪明,等. SAR 图像多时相组合滤波方法[J].计算机工程,2011, 37(21):3-5

[10]Xia Ye. Bam Earthquake: D-INSAR Preliminary Result.[EB/OL]. http://www.gfz-potsdam.de/portal/gfz/Struktur/Departments/Department+1/sec14/topics/radar/sar/earthquake_bam,2012-11-23

[11]Fialko Y, Sandwell D, Simons M, et al. Three-dimensional Deformation Caused by the Bam, Iran, Earthquake and the Origin of Shallow Slip Deficit[J]. Nature,2005, 435(7040): 295-299