合成孔径雷达干涉测量(interferometric synthetic aperture radar, InSAR)是一种大地测量技术，能够全天候获取高精度、连续覆盖的地面高程和地表形变信息。InSAR已在地形测绘、全球环境变化、灾?#981;嗖馄拦赖认喙亓煊虻玫搅斯惴河τ貌⑷〉昧艘幌盗谐晒?sup>^[1]。

PS-InSAR(persistent scatterer InSAR)是一种时间序列InSAR技术，利用多景同一地区的SAR影像，通过统计分析时间序列上幅度和相位信息的稳定性，探测不受时间、空间基线去相关影响的稳定目标。这些目标可能是人工建筑物、裸露的岩石、人工布设的角反射等，由于它们在时间序列SAR影像中几乎不受斑点噪声的影响，经过很长时间间隔仍然保持稳定的散射特性^[2]。

使用PS-InSAR技术能够进行大范围、无接触、点状的mm级地表形变监测。朱叶飞等使用PS-InSAR技术获得了1995-2000年苏州地区地面沉降场^[3]；汪宝存等进行了郑州市地面沉降监测的调查研究^[4]；雷坤超等提取了廊坊市城区2003-2007年时间序列地面沉降形变信息^[5]；张学东等提取了京沪公路(北京-河北)沿线的沉降速率图和沉降剖面图^[6]。

当前，PS-InSAR成果较常见的展示方法是使用专题地图分级设色的方式进行展示，这种方式能够反映PS点的空间分布和形变速率大小两个方面的信息，但是由于二维展示的局限性，不能呈现高程方面的信息，而PS-InSAR数据形成的PS点所包含的空间信息中除了有平面的x、y信息外，还包含高程z的信息，这使得二维地图展示丢失了高程信息维度。

Google Earth是一款强大的三维地图软件，拥有较大覆盖范围的高分辨率影像。段慧娟等引入了海底地形图，并对海底地形和重力异常值进行了标注，为实现海洋物理场数据的科学管理和可视化提供了有效的方法^[7]。王利锋等基于Google Earth对遵化市D级GPS控制点信息进行导入与地标生成，可以实时动态地管理和利用控制点信息，方便了控制点在测绘领域的应用^[8]。黄亚峰等实现了震中分布图、地质资料图、地震现场信息在Google Earth平台上的融合显示^[9]。张荐铭等提出一套基于Google Earth组件的开发方法，设计并实现了天然气加气站信息管理系统，且具备一定的空间分析能力^[10]。

Google Earth可为PS-InSAR成果可视化展示提供良好的场景，能够补充专题地图在高程维度上信息的欠缺，提升展示能力，并且由于其具有更加良好的交互性，用户通过交互操作，能够查看单个PS点的形变速率及累积形变量信息。

本文使用Google Earth作为PS-InSAR监测成果的展示场景，通过构建不同细节层次(level of detail，LoD)KML(keyhole markup language)文件进行展示。

1 PS-InSAR地表形变监测可视化展示方法 1.1 PS-InSAR数据

PS-InSAR是针对高相干PS点的相位处理技术，PS点的差分干涉相位表达式为：

$ {\phi _{{\rm{diff}}}}{\rm{ = }}{\phi _{{\rm{flat\_error}}}}{\rm{ + }}{\phi _{{\rm{topo\_error}}}}{\rm{ + }}{\phi _{{\rm{def}}}}{\rm{ + }}{\phi _{{\rm{atm}}}}{\rm{ + }}{\phi _{{\rm{noi}}}} $

(1)

式中，ϕ_diff为差分干涉相位；ϕ_{flat_error}为卫星轨道数据误差所导致的残余平地相位；ϕ_{topo_error}为外部DEM(digital elevation model)误差所引起的地形相位误差；ϕ_def为两次成像器件的地表形变相位；ϕ_atm为两次大气介质不均一性带来的大气相位；ϕ_noi为噪声相位。

$ {\phi _{{\rm{topo\_error}}}} =-\frac{{4{\rm{ \mathit{ π} }}{B_ \bot }}}{{\lambda R{\rm{sin}}\theta }}\Delta {h_{{\rm{error}}}} $

(2)

$ {\phi _{{\rm{def}}}} =-\frac{{4{\rm{ \mathit{ π} }}{B_ \bot }}}{\lambda }\delta d $

(3)

$ {\phi _{{\rm{diff}}}} =-\frac{{4{\rm{ \mathit{ π} }}{B_ \bot }}}{{\lambda R{\rm{sin}}\theta }}\Delta {h_{{\rm{error}}}} + \frac{{4{\rm{ \mathit{ π} }}}}{\lambda }\delta d + {\phi _{{\rm{res}}}} $

(4)

$ {\phi _{{\rm{res}}}} = {\phi _{{\rm{flat\_error}}}} + {\phi _{{\rm{atm}}}} + {\phi _{{\rm{noi}}}} $

(5)

假设地表形变为线性形变，则有δd=vT，v为形变速率，综合上述公式，有：

$ {\phi _{{\rm{diff}}}} = {k_1}{B_ \bot }\Delta {h_{{\rm{error}}}} + {k_2}vT + {\phi _{{\rm{res}}}} $

(6)

$ {k_1} = \frac{{4{\rm{ \mathit{ π} }}}}{{\lambda R{\rm{sin}}\theta }}, {k_2} =-\frac{{4{\rm{ \mathit{ π} }}}}{\lambda } $

(7)

$ {\phi _{{\rm{diff}}}}^{i-{\rm{ref}}} = \Delta {k_1}^{i-{\rm{ref}}}{B_ \bot }\delta \Delta {h_{{\rm{error}}}}^{i-{\rm{ref}}} + \Delta {k_2}^{i - {\rm{ref}}}\delta \Delta {v_{{\rm{error}}}}^{i - {\rm{ref}}}T + {\phi _{{\rm{pres}}}}^{i - {\rm{ref}}} $

(8)

式中，ϕ_diff^i-ref、δΔh_error^i-ref、δΔv^i-ref分别为PS点对之间的差分相位之差、DEM误差校正之差和形变速率之差。通过回归分析可得到形变速率和高程改正值，高程改正值叠加DEM信息为PS点高程信息，高程信息和平面坐标信息构成了PS点的三维空间信息。PS-InSAR需要展示的信息如表 1所示。

1.2 KML数据构建

KML是一种开放标准，其官方名称为OpenGIS KML编码标准(OpenGIS KML encoding standard，OGC KML)，KML用于在地球浏览器(如Google地球、Google地图和Google地图移动版)中显示地理数据^[11]。其使用包含嵌套的元素和属性的结构(基于标记)，并符合XML标准。由于PS-InSAR监测结果通常会产生大量的数据点，在较高的浏览视角进行全局浏览时，会同时载入过多的点，从而造成软件响应缓慢。因此，本文使用KML地图项“区域”进行展示数据的构建，只有当用户位于视图内且占用一定的屏幕部分时，区域才会加载并绘制该数据。因此，可在Google Earth中添加大型数据集，而不会降低其性能。

区域的构建主要包含如下属性设置：边框、细节层次级别、海拔、嵌套区域、基于区域的NetworkLink等。其中，边框<LatLonAltBox>是包括一组对象或数据点的范围，包含东、西、南、北边界以及最低海拔<minAltitude>和最高海拔<maxAltitude>，如图 1所示。

在区域中，<minLodPixels>和<maxLodPixels>元素可指定屏幕的区域(以方形像素为单位)。当数据投影到屏幕上时，所占用的屏幕区域必须大于<minLodPixels>且小于<maxLodPixels>，数据才可显示。区域的投影大小超出这些限制后，投影将不再可见，区域会变为非激活状态。

本文通过分级方式进行KML文件建立，层级由低到高，PS点数增加。最高级的KML包含全部PS点，对其他级别的PS数据进行点抽稀，抽稀时可对形变速率进行排序，对形变速率较大的点予以保留，突出沉降区域的细节。同时，对每一级利用“区域”分块构建KML数据，级数越高，分块越多。分级区块构建完成后，使用NetworkLink将各个区块的KML连接起来，得到最终展示数据。构建KML文件流程如图 2所示。

2 实验及分析

本文使用四川省地理省情监测成都市地表形变监测数据作为实验数据进行构建，监测时间为2013-04~2015-01。选取天府广场周边共计34 518个PS点。所用PS-InSAR成果数据构成如图 3所示，其中包含形变速率(mm/a)、空间x、y、z坐标(经度、纬度、高程)、各期累积形变量(mm)等3个方面的信息，数据使用二维表格进行存储。

KML文件符合XML标准，针对KML文件，可按照操作XML的方式进行。因此，本文使用C#语言编写程序，构建各级KML数据。初始KML文件可使用Google Earth或GMT Tools创建，然后根据需要增删所用标签及内容。

在构建的3个级别的KML文件中，第三级不经筛选，选用全部PS点，最高级别各区域中PS点数不宜过多(不超过2 000个)。一、二级逐级抽稀后，再细分为多个区域进行构建，数据构建完成后加载进入Google Earth中进行展示，优化“区域”所涉及的各个标签参数设定。各级区域示意图见图 4。

所构建的KML文件主要分为节点文件、数据文件、样式文件、图例文件等4种功能类型，节点文件(一级、二级、三级)不含数据内容，使用<Networklink>进行文件的逐级连接，数据节点文件中包含点位的坐标和累积形变量信息。4种文件的标签构成如图 5所示。其中，图 5(a)为节点KML文件，主要包含NetworkLink，通过<href>标签连接一级、二级、三级KML文件；图 5(b)为数据文件，主要包含Placemark标签，空间坐标写入Placemark>Point>coordinates标签，形变速率数据和累积形变量数据写入Placemark>ExtendedData>Data>value标签；图 5(c)为Style文件，采用分级设色进行形变速率颜色的渲染，共构建15个Style，用ID属性和Placemark中的styleUrl连接，使用JavaScript语言在Style>BalloonStyle>text标签中编写脚本，用于构建折线图，展示累积形变量信息; 图 5(d)为Legend标签用overlayXY、screenXY设置图例显示位置，使用Icon>href标签设置图例(位图文件)的硬盘位置。

构建完成的KML文件结构如图 6所示。

在完成构建后，将CDMain.kml文件加载入Google Earth中进行展示，如图 7所示。其中，图 7(a)为全局浏览时显示的第一级概略数据；图 7(b)为更为精细的数据显示；图 7(c)为天府广场东俯瞰图；图 7(d)为四川省科技馆形变情况俯瞰图；图 7(e)为城市之心大厦周边形变俯瞰图；图 7(f)为单击形变点后展示的该点累积形变量折线图。

从图 7可以看出，利用Google Earth作为展示场景，构建KML文件进行PS-InSAR形变结果展示，相对于专题地图进行展示的方式，能够展示PS点的三维空间信息，具有良好的交互性。

3 结束语

本文利用Google Earth作为PS-InSAR形变监测成果展示场景，使用示例PS-InSAR数据构建了KML文件，并进行了展示。结果表明，使用KML文件在Google Earth中进行展示，能够在3个空间维度上进行PS点的绘制，解决了专题地图将所有PS点投影在二维平面上造成的信息缺失问题。