引用本文
彭鹏, 杨红磊. 基于改进相干点目标技术的阜阳市地面沉降调查[J].中国地质调查, 2016,3(4): 63-68.
PENG Peng, YANG Honglei. Study of land subsidence monitoring based on an improved coherent point target technology in Fuyang City[J].Geological Survey of China, 2016,3(4): 63-68.
  
Permissions
Copyright@2016, 《中国地质调查》编辑部
版权所有 @ 《中国地质调查》编辑部
基于改进相干点目标技术的阜阳市地面沉降调查
彭鹏1, 杨红磊2
1.安徽省地质调查院,合肥 230001
2.中国地质大学(北京),北京 100083

第一作者简介: 彭鹏(1982—),男,工程师,硕士研究生,主要从事遥感地质调查工作。Email: 119170256@qq.com。

收稿日期: 2016-03-10
基金资助: 中国地质调查局 “全国地表形变遥感地质调查(安徽省)(编号: 1212011140037)”和安徽省国土资源厅“基于PS-InSAR的皖北地区长时序地面沉降监测方法研究(编号: 2016-k-8)”项目共同资助
摘要

提出了一种改进的相干点目标时序InSAR技术方法。该方法基于子视相干值识别小数据集的相干点目标,采用趋势项滤波方法去除大尺度残差相位,并以阜阳市2012—2013年间10景Radarsat-2 Wide模式数据为例,获得了高空间密度的时序高相干点目标,提高了结果的可靠性,获取的沉降区数据与前人实测资料及野外调查的沉降现象相一致,为开展阜阳市地面沉降风险分析和经常性预警监测提供技术方法和示范应用。

关键词: 相干点目标技术; 地面沉降; 滤波; Radarsat-2
中图分类号:TP79 文献标志码:A 文章编号:2095-8706(2016)04-0063-06
Study of land subsidence monitoring based on an improved coherent point target technology in Fuyang City
PENG Peng1, YANG Honglei2
1. Geological Survey Institution of Anhui Province, Hefei 230001, China
2. China University of Geosciences(Beijing), Beijing 100083, China
Abstract

This paper proposed an improved coherent point target InSAR technology. The value of the sub-look coherent was used to identify the coherent point targets of small dataset, and the trend filter was then applied to remove the large-scale residual phase. Using 10 Radarsat-2 Wide Model images over Fuyang City during 2012—2013, the high density coherent point targets were obtained and the reliability of the results were improved. The obtained land subsidence was in accordance with the previous measurement data and our field survey. This method provides the technical supports and application example for risk analysis and the regular monitoring of the land subsidence in Fuyang City.

Keyword: coherent point target technology; land subsidence; filter; Radarsat-2
0 引言

InSAR(Interferometric SAR)技术是利用不同时期2景SAR影像间的相位差来获取地面沉降信息, 具有全天时、全天候和高形变监测精度的特点[1]。利用该技术已在我国主要地表形变区和重大工程区开展了大量地表形变摸底调查工作[2]。相干点目标技术(coherent point target, CPT)是通过对时间序列中相位保持稳定的点进行建模分析, 获得高精度的地表形变信息, 被广泛应用在地面沉降[3, 4]、地震形变[5]、火山运动[6]以及山体滑坡[7]等领域。目前识别高质量相干点目标的常用方法有相干系数阈值法[8]、振幅离差阈值法[8]、相位离差阈值法[9]和信噪比[10]等。然而, 上述方法都基于大数据量(> 30景), 从统计的思想选择相干点目标, 不适合小数据集的情况。虽然相干点目标分析(interferometric point target analysis, IPTA)技术采用了子视相干值来识别小数据集的相干点目标[3], 但在稳健性方面仍有不足。针对这一问题, 提出了利用子视相干值离差识别相干点目标的方法。数据源采用Radarsat-2数据, Radarsat-2 Wide模式具有监测范围大的特点, 适用于地面沉降普查和经常性日常监测预警。但是由于其定轨精度较低, 采用常规的FFT[11]、GCP[12]和二项式拟合方法[13]也难以消除轨道残余误差, 基于该误差在空间上表现为大尺度低频信号[14], 提出了采用趋势项滤波方法纠正其影响。

1 研究方法

改进的相干点目标技术共包含识别相干点目标、生成差分干涉图、相位解缠和估计形变值4个部分。

1.1 识别相干点目标

首先从若干景SAR影像中, 依据总体相关性最大原则[14]确定参考影像, 然后把其他影像配准到该影像的成像几何, 配准精度优于0.2个像元。

与分布式散射体(如农田、森林等)相比, 点状目标(如建筑物、船舶等)可视为由一个或者多个“ 主散射体” 构成的反射体, 具有强振幅值, 不受斑点噪声的影响, 在单视复数影像中具有确定的相位值, 在频率域中表现为一种独特的光谱特征, 可以依据该属性来识别相干点目标。通过在方位向或者距离向降低处理器带宽, 将频谱分割为若干部分, 即为子视分解[15]。由于分布式散射体内各独立散射体回波的随机性, 当任意子视影像的频谱不重叠时, 其斑点噪声像素之间是失相干的, 可以采用相关系数γ [16]确定相干散射体, 即

γ = < X1X2* > |< X1X1* > < X2X2* > , 0≤ γ ≤ 1 。 (1)

式中: X1、X2分别为子视光谱; * 代表取共轭; < > 代表空间求平均。依据式(1)采用单景SAR影像就可以有效地确定相干点目标, 但是其中包含了汽车、轮船等移动目标。在CPT处理中, 仅需要在整个时间序列内相对稳定的相干点目标, 故提出了采用光谱离差DC衡量时间序列点的光谱稳定性, 剔除移动相干点目标, 即

DC= δCmC, (2)

式中δ C和mC分别为时间序列内光谱相关值的标准差和均值。

1.2 生成差分干涉图

为了弱化时空失相干的影响, 采用短时空基线的原则确定干涉对。对于不同的波长, 设置不同的基线阈值, 一般建议设为极限基线的1/3。干涉组合确定以后, 对相干点目标进行干涉, 得到每个点的干涉相位。再依据轨道和地形数据模拟平地效应和地形相位, 并从干涉图中去除, 得到每个点的差分干涉相位。此时的差分干涉相位包含轨道残余相位、地形残余相位、形变相位、大气延迟相位和噪声相位。除了地形残余相位和噪声相位外, 其他部分都具有空间相关性。

1.3 相位解缠

由1.2节获得的相干点目标x的一次差分时间序列干涉相位为

Φxk= φx, topok+ φx, defok+ φx, orbitk+ φx, atmok+ φx, noisek, (3)

式中: k为干涉图序列号; φx, topok为地形残余引起的相位值; φx, defok为地物在雷达视线方向移动引起的相位值; φx, orbitk为轨道残余误差相位; φx, atmok为大气延迟引起的相位; φx, noisek为失相关噪声相位。

从相干点目标集中选择一个具有稳定高相干值且位于地表比较稳定区域的点作为计算参考点y, x和y再次做差分, 得到二次差分干涉相位, 表达式为

Φx, yk= φ(x, y), topok+ φ(x, y), defok+ φ(x, y), atmok+ φ(x, y), orbitk+ φ(x, y), noisek, (4)

式中: φ(x, y), topok为相对参考点高程残余差引起的地形残余相位值; φ(x, y), defok为相对参考点视线方向位移差对应的相位值; φ(x, y), atmok为相对参考点的大气延迟相位差; φ(x, y), orbitk为相对参考点的轨道残余误差相位; φ(x, y), noisek为相对参考点的失相关噪声相位差。

由于相邻点之间的高程差是垂直基线的函数, 相邻点之间的形变速率与时间间隔相关, 式(4)也可以表示为

Φx, yk=[β kΔ h(x, y)+ 4πλTkΔ v(x, y)]+ φ(x, y), atmok+ φ(x, y), orbitk+ φ(x, y), noisek, (5)

式中: λ 为雷达波长; β k为高程常数; Δ h(x, y)为相邻点高程差; Tk为相邻时刻SAR时间间隔; Δ v(x, y)为时序形变速率。

采用解空间搜索法求得高程差与形变速率, 恢复相邻点之间在每个时序差分干涉对上的真实相位, 实现相干点目标时间域相位解缠。从原始差分干涉相位中减去高程和形变相位后, 得到的残差相位包含轨道残余误差、大气延迟、非线性形变和噪声相位。由于轨道残余误差相位和大气延迟相位在空间上具有强相关性, 对解缠后的残差相位采用趋势项滤波方法进行分离, 该方法由多视平均、空间滤波和空间插值3部分组成。具体流程如下:

(1) 多视平均。假设M× N的矩阵D对应的矩阵元素为di, j, 采用M× N的多视窗口, 得到多视平均后的结果为

d—= ∑i=1M∑j=1Ndi, jMN。 (6)

由式(6)可知, 对数据进行多视平均处理, 去除噪声的同时, 还保持了趋势项相位的形态。

(2) 空间滤波。由于非线性形变相位和轨道残余相位、大气延迟相位一样, 在空间上都具有较强相关性, 但是空间尺度较小, 当采用视数较大的窗口多视处理后, 非线性形变表现为噪声点, 可以采用空间均值、空间中值等滤波方法进行去除。

(3) 空间插值。空间滤波后轨道残余相位和大气延迟相位的空间分辨率和原始图像不一致, 需要采用空间插值方法(反距离加权、样条函数插值和Kriging插值等)恢复到原始分辨率。

1.4 估计形变值

1.3节的处理不仅实现了差分干涉相位的相位解缠, 而且可以分离出轨道残余相位和大气延迟相位。从得到的解缠差分干涉图中剔除轨道残余相位和大气延迟相位后, 剩下的差分干涉相位只包含了线性形变、非线性形变和高程残差信息, 采用SVD[17]分解便可以获得时间序列形变的最小二乘解。

2 阜阳市地面沉降监测应用

阜阳市位于安徽省西北部, 黄淮海冲积平原南端, 地形平坦开阔。境内西北高、东南低, 高程为2039 m, 坡降为1/8 0001/9 000。利用InSAR数据尚未开展过系统的地面沉降调查工作。根据阜阳市一、二等水准测量结果插值绘制的1998— 2008年10 a间的沉降等值线图, 沉降中心位于三角游园一带, 沉降速率可达-40 mm/a[18], 近些年由于限制地下水开采, 沉降变化缓慢。

数据源为Radarsat-2 Wide模式数据, 空间分辨率距离向为11.8 m, 方位向为5.2 m, 极化方式为VV, 获取时间为2012年2月4日— 2013年1月29日, 共10景。高程数据采用SRTM DEM。

依照总体相干值选择2012年5月10日为参考影像, 把其他影像配准到该影像的成像几何。分别对每景SAR影像进行子视分解, 求取子视相干值, 设置时序子视相干均值阈值为0.37, 相干离差阈值为1.15, 共识别 673 456 个相干点目标, 主要

分布在阜阳市区及郊区的居民地, 农田区域也有零星分布, 实地调查得知, 主要是农田中的灌溉设施。设置空间基线阈值为300 m, 时间基线为200 d, 共确定30对干涉对, 干涉对组合如图1所示。

图1 干涉对基线组合示意图Fig.1 Schematic diagram of interferometer baseline combination

由2012年2月4日与2012年4月16日数据组合干涉对生成的差分干涉图如图2所示。

图2 2012-02-04与2012-04-16影像的差分干涉图Fig.2 Differential interferograms of images on 2012-02-04 and 2012-04-16

从图2(a)可以看出, 由于轨道定轨精度低, 生成的基线误差较大, 造成初始差分干涉图中残余的平地效应和高程误差较大, 甚至淹没了形变值对应的相位。采用GCP纠正虽然能够消除一部分基线误差(图2(b)), 但是仍不能完全消除其影响。图2(c)为采用趋势项滤波方法得到的改进后差分干涉图(多视视数为10× 10, 空间滤波方法为均值滤波, 滤波窗口为7× 7, 插值方法为Kriging插值), 可以明显发现已基本消除大尺度残余相位的影响, 能够识别出形变区域, 说明本文方法是可行的。

对所有干涉对成功相位解缠后, 采用SVD分解得到沉降速率如图3所示。

由图3可知, 阜阳市形变主要集中在阜阳市区, 存在一个较明显的沉降滤斗, 市区南部村镇也有沉降现象。该漏斗以奎星楼为中心, 东至三角洲公园、西至中南现代城、南至区政府、北至青颖公园, 覆盖范围约14 km2, 最大形变速率为21 mm/a, 分布范围与以往沉降中心一致。由于2008年后该地区限制开采地下水, 沉降有所减缓, 本次InSAR调查速率比参考文献[18]介绍的有所减小是合理的。

图3 阜阳市2012— 2013年间地面沉降速率Fig.3 Land subsidence rates in Fuyang city from 2012 to 2013

3 结果验证与机理分析

通过分析和野外调查发现, 阜阳市和地面沉降直接相关的现象有抽水井管道弯曲、抬升和桥梁拉伸、错位。

3.1 井管弯曲、抬升地表形变现象

开采地下水的深井井管抬升及井台地面开裂等地表变形现象如图4和图5所示。

图4 自来水公司深井井管弯曲Fig.4 Bended pipe of bore in the Fuyang tap water company

图5 井台遭受破坏Fig.5 Damaged well platform

具体统计数据如表1。

表1 由地面沉降引发变形井的统计数据 Tab.1 Statistical data of bended well caused by land subsidence

由图4、图5和表1可以看到, 位于InSAR地面沉降区的服装厂深井(现属自来水公司), 井台较周围地面有明显抬升, 井台地面开裂, 裂缝呈放射状; 该深井井管弯曲明显, 可达25 cm。该形变主要由于地下水超采造成浅部地层沉降所致, 而井管末端所在深部地层较为稳定, 使最上端的井台与井管保持在原始高度, 产生类似井管蹿升的效果。野外调查点均位于InSAR形变速率图中存在明显沉降的区域。

3.2 桥梁拉伸、错位地表变形现象

阜阳市跨度较大的大桥、公路等都出现了桥面开裂、拉伸及错位等现象(如图6、7所示), 经过修缮, 已恢复使用。

图6 桥梁错位Fig.6 Bridge dislocation

图7 桥梁拉伸Fig.7 Bridge stretching

通过对InSAR沉降区进行野外实地调查, 该现象都位于沉降区范围, 可能为不均匀沉降所致, 应作为地面沉降重点调查区。

3.3 地面沉降形成机理分析

3.3.1 地层特征

研究区区域地层属华北地层区淮河地层分区阜阳地层小区, 发育有新元古界霍邱群和五河群, 古生界寒武系、奥陶系、石炭系和二叠系, 中生界白垩系, 新生界古近系、新近系和第四系地层。研究区内主要为新生代松散层覆盖, 厚800900 m, 其中第四系厚130150 m, 松散沉积巨厚。

3.3.2 水文地质特征分析

阜阳市地下水类型为单一松散岩类孔隙水, 根据地下水赋存介质的性质、类别和组合的不同, 自上而下划分为浅层地下水和深层地下水[18, 19] 。浅层地下水赋存于50 m以浅的全新世、晚更新世地层中, 与大气降水、地表水关系密切, 按上下关系可称其为第一含水层组(浅层); 深层地下水赋存于50 m以下的地层中, 与大气降水、地表水关系不密切。根据水文地质结构和目前开采现状, 将深层地下水进一步划分为2个含水层组, 即第二含水层组(中深层, 埋深50150 m)和第三含水层组(深层, 埋深150500 m)。随着城市化的发展, 深井越打越多, 越打越深, 开采量不断增加, 致使城区及周边深层地下水超采, 这也成为研究区发生地面沉降的主要原因。

阜阳市含水层模型剖面示意图如图8所示。

图8 阜阳市含水层模型剖面示意图Fig.8 Schematic section of aquifer model in Fuyang

4 结论

本文分别采用子视相干点目标识别方法和趋势项滤波方法, 对传统相干点目标技术中的相干点目标识别和大尺度残余误差剔除2方面进行改进。以阜阳市2012— 2013年间10景Radarsat-2 Wide模式数据为例进行实验, 获得该市2012— 2013年间的形变速率图, 认为该市存在一个较为明显的形变漏斗。结合形变区实地调查, 初步分析了形变产生的影响及原因, 为该市形变监测网的布设和地面沉降防治提供了决策依据。

(责任编辑: 刁淑娟)

The authors have declared that no competing interests exist.

参考文献
[1] Gabriel A K, Goldstein R M. Crossed orbit interferometry: Theory and experimental results from SIR-B[J]. International Journal of Remote Sensing, 1988, 9(8): 857-872. [本文引用:1]
[2] 范景辉, 燕云鹏, 葛大庆, 等. 全国地表形变遥感地质(InSAR)调查技术指南[M]. 北京: 地质出版社, 2015. [本文引用:1]
[3] Werner C, Wegmüller U, Strozzi T, et al. Interferometric point target analysis for deformation mapping[C]//Proceedings of 2003 IEEE International Geoscience and Remote Sensing Symposium. Toulouse, France: IEEE, 2003, 7: 4362-4364. [本文引用:2]
[4] Yang H L, Peng J H. Monitoring urban subsidence with multi-master Radar interferometry based on coherent targets[J]. Journal of the Indian Society of Remote Sensing, 2015, 43(3): 529-538. [本文引用:1]
[5] Feng G C. Coseismic deformation and ionospheric variation associated with Wenchuan earthquake estimated from InSAR[D]. Hong Kong: The Hong Kong Polytechnic University, 2011. [本文引用:1]
[6] Lu Z, Mann D, Freymueller J T, et al. Synthetic aperture Radar interferometry of Okmok volcano, Alaska: Radar observations[J]. Journal of Geophysical Research, 2000, 105(B5): 10791-10806. [本文引用:1]
[7] Berardino P, Costantini M, Franceschetti G, et al. Use of differential SAR interferometry in monitoring and modelling large slope instability at Maratea (Basilicata, Italy)[J]. Engineering Geology, 2003, 68(1/2): 31-51. [本文引用:1]
[8] Ferretti A, Prati C, Rocca F. Permanent scatterers in SAR interferometry[J]. IEEE Transactions on Geoscience and Remote Sensing, 2001, 39(1): 8-20. [本文引用:2]
[9] Ferretti A, Prati C, Rocca F. Nonlinear subsidence rate estimation using permanent scatterers in differential SAR interferometry[J]. IEEE Transactions on Geoscience and Remote Sensing, 2000, 38(5): 2202-2212. [本文引用:1]
[10] Adam N, Kampes B M, Eineder M. Development of a scientific permanent scatterer system: Modifications for mixed ERS/Envisat time series[C]//Proceedings of the 2004 ENVISAT and ERS Symposium. Salzburg: ADS, 2004. [本文引用:1]
[11] Hanssen R F. Radar Interferometry: Data Interpretation and Error Analysis[M]. Netherland s: Springer, 2001. [本文引用:1]
[12] 陈刚, 汤晓涛, 余旭初, 等. 基于GCP的星载InSAR系统基线估计算法精度分析[J]. 测绘通报, 2009(6): 16-18, 21. [本文引用:1]
[13] 杨红磊, 彭军还, 张丁轩, 等. 轨道误差对InSAR数据处理的影响[J]. 测绘科学技术学报, 2012, 29(2): 118-121, 126. [本文引用:1]
[14] Zebker H A, Villasenor J. Decorrelation in interferometric Radar echoes[J]. IEEE Transactions on Geoscience and Remote Sensing, 1992, 30(5): 950-959. [本文引用:2]
[15] Cumming I G, Wong F H. 合成孔径雷达成像——算法与实现[M]. 洪文, 胡东辉, 译. 北京: 电子技术出版社, 2007. [本文引用:1]
[16] Souyris J C, Henry C, Adragna F. On the use of complex SAR image spectral analysis for target detection: Assessment of polarimetry[J]. IEEE Transactions on Geoscience and Remote Sensing, 2003, 41(12): 2725-2734. [本文引用:1]
[17] Berardino P, Fornaro G, Lanari R, et al. A new algorithm for surface deformation monitoring based on small baseline differential SAR interferograms[J]. IEEE Transactions on Geoscience and Remote Sensing, 2002, 40(11): 2375-2383. [本文引用:1]
[18] 安徽省地质调查院. 安徽省阜阳市地面沉降调查报告[R]. 2004. [本文引用:2]
[19] 安徽省地矿局第一水文地质工程地质队. 安徽省阜阳市水、工、环地质综合详查报告(1/5万)[R]. 1988—1990. [本文引用:1]