由于对石油、天然气、地下水、卤水等自然资源的长期大量开采,以及地表荷载增加、新近沉积物固结压实等因素的影响,现代黄河三角洲地区地面沉降现象已经引起广泛关注^[1-8]。地面沉降导致建筑物、公共基础设施、输油管线破坏,海平面相对上升,海水入侵和风暴潮灾害加剧^[2,3,5,7]。因此,对现代黄河三角洲地面沉降进行监测显得越来越重要。由于早期黄河三角洲地面沉降监测站少,覆盖面不全,监测网密度较低,很难准确地定义当时的地面沉降区的范围和形变幅度,这势必影响对黄河三角洲地面沉降情况进行全面系统的认识和研究,因此需要开展监测工作。

合成孔径雷达干涉测量（Interferometric Synthetic Aperture Radar,InSAR）技术是利用对同一地区观测的两景或多景SAR影像复数据进行相干处理,获取高精度的三维地形信息及微小形变信息的技术^[9-11]。InSAR技术能够全天候、全天时、大覆盖范围、低成本、高时间分辨率、高空间分辨率、高精度（理论上能达到mm级）的完成对地表的形变监测。2001年A Ferretti等人提出PSInSAR（Permanent Scatters InSAR）技术^[12-14],又称永久散射体技术。永久散射体是指在相当长时间跨度InSAR影像序列中保持稳定且较高后向散射特性的目标点,他们具有较高相干性,受空间、时间影响比较小。PSInSAR技术监测缓慢地表形变具有无与伦比的优势,为地面沉降的预防和治理提供重要的基础数据。本文针对黄河三角洲地区湿地广泛分布和植被茂密的特点展开研究,对PSInSAR技术在复杂环境中的应用进行实验。此外,由于传统的PSInSAR技术对数据量的要求较高,一定程度上限制了PSInSAR技术的应用,本文利用1992年6月15日至1995年8月2日期间覆盖现代黄河三角洲部分地区的10景ERS1/2影像开展基于小数据集的PSInSAR技术监测地面沉降应用研究,并利用水准观测数据对PSInSAR技术获得的结果进行验证。

1 研究区概况

现代黄河三角洲是1855年以来,黄河冲淤作用形成的冲积扇,是最年轻的陆地。研究区位于莱州湾的西岸（图1）,包括垦利县、东营市、东营区及周边乡镇等。该区域地形平坦,海拔高度多小于6 m,平均坡降平缓,为1/10 000左右。研究区范围为118°~119°E,37°~38°N,覆盖面积约70 km×70 km。研究区内石油、天然气、卤水和地热资源等资源丰富。研究区植被茂密,湿地广泛分布,这使得传统测量中观测点的设置、维护和观测相对困难。

显示原图| 下载原图ZIP| 生成PPT

图 1   研究区地理位置

Fig.1   The location of study area

2 PSInSAR技术原理

Hooper提出基于相位特征的PS算法的基本原理如下^[14]：获取覆盖研究区的N+1景SAR影像。首先将辅影像与主影像进行配准,干涉处理得到N景差分干涉图,利用外部DEM去除地形相位。采用基于幅度分布指数 DA（Amplitude Dispersion Index）来预先选定候选PS。得到PS候选点后,分析其相位稳定性筛选最佳的PS点,基本模型如下:

第i景干涉图中第x个PS点的相位残差是：

ϕx,i=ϕdef,x,i+ϕa,x,i+ϕorb,x,i+ϕε,x,i+nx,i（1）

其中, ϕx,i是第i景干涉图的第x个像素的相位值, ϕdef,x,i是雷达视线方向的形变相位, ϕa,x,i是卫星不同时间过境大气延迟相位差, ϕorb,x,i是轨道误差相位, ϕε,x,i是DEM误差相位（与垂直基线B_n-i成比例）, nx,i是噪声相位,包括像素内散射体变化、热效应、配准误差等引入的噪声相位。

假设 ϕdef, ϕα和 ϕorb在一定长度L范围内是空间相关的, ϕε和 n在这个范围内是不相关的,均值为零。以像素X为圆心,以L为半径的一个圆形区域内的所有相位的均值可表示为：

ϕ̅x,i=ϕ̅def,x,i+ϕ̅α,x,i+ϕ̅orb,x,i（2）

公式（1）减去（2）得到：

ϕx,i-ϕ̅x,i=B⊥,x,iKε,x+nx,i-n̅'x,i（3）

公式（3）中, n̅'=(ϕ̅def-ϕdef)+(ϕ̅α-ϕα)+(ϕ̅orb-ϕorb)。DEM误差相位 ϕε,x,i和垂直基线 B⊥,x,i成比例关系, Kε,x为比例系数。

利用所有可以利用的干涉图,对像素x使用最小二乘法来估计比例系数 Kε。基于像素x的时间相干性 γx定义为：

γx=1N∑i=1Nexp{j(ϕx,i-ϕ̅x,i-ϕ^ε,x,i)}（4）

其中,N是使用的干涉图的数量, ϕ^ε,x,i是 ϕε,x,i的估计值,假设 n̅'ε,x,i很小, γx是像素的相位稳定量,也是判断PS点的指标。

3 数据处理

共收集了1992年6月15日至1995年5月24日期间覆盖现代黄河三角洲地区的10景ERS1/2降轨影像进行实验。综合考虑了时间基线、空间垂直基线和多普勒质心频率基线三种基线的临界值和三种基线绝对值之和最小两种因素,选择1993年11月22日的影像为主影像。去除地形相位时采用外部DEM,选用覆盖研究区域的SRTM数据,该数据平面分辨率为90 m,在地形起伏较小的地区高程精度较高。

利用Hooper提出的新的PS算法对数据进行处理,主要步骤如图2所示：首先,选取主影像,裁剪研究区,并将主影像与辅影像依次配准并干涉处理,去除地形相位影响,得到9景差分干涉相位图;其次,确定PS点,将幅度分布指数设定为阈值0.4筛选PS候选点,并对空间相关误差和空间非相关视角误差进行估计,通过估计PS候选点的相位稳定性确定最终的PS点;然后去除估计的误差并利用三维相位解缠方法进行相位解缠;最后,在时间域上做高通滤波,在空间域上做低通滤波,去除残余DEM误差、大气干扰与轨道误差影响,获得研究区的形变相位。

显示原图| 下载原图ZIP| 生成PPT

图2   PSInSAR技术数据处理流程

Fig. 2   Processing flow chart of PSInSAR technical data

4 结果分析

通过对覆盖现代黄河三角洲部分地区的ERS1/2数据进行处理取得了较好结果,有8×10⁴多个相位稳定PS点被选取出来,密度约为20个/km²。为提高监测结果的可靠性,我们以3.0 mm/a为阈值将标准差较大的PS点剔除,最后获得了约46 000个可靠的PS点的形变数据。由于东营机场（见图3中矩形区域）于1985年10月建成,在建设过程中对地面进行了夯实处理,并且建成后至2001年前没有投入使用,因此,不考虑大地构造运动的作用,以东营机场作为参考点,对整个研究区的监测结果进行校正,获得现代黄河三角洲1992年6月15日至1995年5月24日间的平均地面沉降沉降速率（文中提到的PSInSAR技术获得的地面沉降速率不作特殊说明都是指的该值）,并与1995年9月18日TM影像第5波段数据进行叠加显示（图3）其中负值表示沉降。

显示原图| 下载原图ZIP| 生成PPT

图 3   沉降速率（mm/a）

Fig. 3   The map of subsided velocities (mm/a)

在图3中,PS点在居民区、道路、水库堤坝等处分布较密,且明显沿道路和堤坝呈线状分布,证明实验选取的PS点的稳定性和可靠性。此外,由于在植被茂密地区和东部湿地分布区存在人工建筑、油井和裸地等,也都顺利提取到PS点。研究区内地面沉降广泛分布（见图3）,重点选取A、B、C三个沉降量相对较大的椭圆形区域进行讨论,其中A区域为胜利采油厂及其周边地区,提取到的2 439个PS点的地面沉降速率的平均值为10.1 mm/a,区域内有大量的PS点的沉降速率在（6~15）mm/a,且有较多PS点的沉降速率大于15 mm/a,最大沉降速率达到22.0 mm/a。B区域位于东营市西城区, 提取到的5 417个PS点的沉降速率的平均值为7.1 mm/a,B区域内有较多PS点的沉降速率分布在9 mm/a左右,只有较少PS点的沉降速率大于12 mm/a,其中,耿井水库南部油井附近监测到的9个PS点的沉降速率的平均值为14.7 mm/a,最大沉降速率达到16.3 mm/a。C是东营市东城区,提取到的1 156个PS点的沉降速率的平均值为7.4 mm/a,其中有较多PS点的沉降速率分布在12 mm/a左右,只有很少PS的沉降速率大于15 mm/a,最大沉降速率达到17.7 mm/a。此外,东营市垦利县郝家镇胜利油田稠油厂附近油库分布广泛,油库附近各PS点的沉降速率的平均值为6.5 mm/a;胜坨镇最大沉降速率为12.8 mm/a;卧铺乡最大沉降速率为17.2 mm/a;黄河口镇最大沉降速率达到15.6 mm/a;羊口镇清河采油厂附近监测到最大沉降速率达到11.5 mm/a;羊口镇地下卤水储量丰富,盐场广泛分布,盐场中提取到的几个PS点中沉降速率最大为16.8 mm/a。另外,沿海地区地面沉降普遍存在,广北水库附近及其以东地区平均地面沉降速率约为4.4 mm/a。

将A、B、C三个沉降比较集中的区域与该地区1992年的居民区分布图、1980年油田分布图与1995年的油田分布图进行叠加分析（图4）。对图4进行统计分析,发现位于1980年与1995年油田分布区和居民区的A区域,平均沉降速率为16.3 mm/a,最大沉降速率达到22.1 mm/a。位于1980年与1995年油田分布区和居民区均分布的C区域,平均沉降速率为12.0 mm/a,最大沉降速率达到12.4 mm/a。位于居民区与1995年油田分布区内的B区域平均沉降速率为7.0 mm/a。位于A区域上方的垦利县城,仅在在地表荷载的作用下的平均沉降速率是6.8 mm/a。位于C区域东北方向的胜利油田东辛采油厂采油五矿处主要是1980年与1995年油田分布区,建筑物较少,在油田开采活动的影响下平均沉降速率为6.6 mm/a。图3中机场东部沿海的广北水库地区仅在自然状态下平均沉降速率为4.4 mm/a。综上可知,由现代黄河三角洲新成陆和冲淤特性决定的,沉降速率与地表载荷增加、石油开采活动和沉积物的自然固结压实具有较强的相关性。

显示原图| 下载原图ZIP| 生成PPT

图 4   居民区、油田和沉降速率分布

Fig. 4   The distribution of residential area, oilfield and subsided velocity

在A、B、C三个区域,分别选取一个能够反映该区域形变特征的沉降量比较大的PS点,依次标注字母A、B、C,并给出各个时间段内PS点的沉降量情况。其中所有沉降量都是相对1993年11月22日获取的主图像的垂直方向的形变量,正数表示该时间段的地面相对于主影像上升,负数表示该时间段地面相对主影像下沉（单位mm）。图5中依次是A、B、C三个区域中选取的PS点的时间序列沉降情况。其中,A区域PS点的沉降量范围是13.3~-27.7 mm;B区域PS点的沉降量范围是10.0~-13.6 mm;C区域PS点的沉降量范围是13.2~-17.4 mm。

显示原图| 下载原图ZIP| 生成PPT

图 5   A、B、C三个PS点形变时间序列

Fig. 5   Time series of A、B、C PS deformation points

利用收集到的4个水准点数据对PSInSAR监测结果进行验证,情况见表1（图3中随着经度的增大依次为水准点1~4）。水准点3在1990~1993年间没有沉降,1993~1994年沉降速率为12 mm/a,提取到的距离最近的PS点相距约80 m,沉降速率为14.5 mm/a,由于相邻时间段沉降速率变化不大,假设1995年初仍然保持12 mm/a的速率匀速下沉,将获得水准观测结果进行插值与PSInSAR技术监测结果进行对比（如图6）,计算获得PSInSAR技术监测结果和水准监测结果间的中误差为0.84 mm。通过以上分析,考虑到位置和观测时间的差异,结合该区域形变缓慢、范围大和相邻时间段沉降速率变化不大的特点,认为PSInSAR技术的监测结果与水准观测结果吻合度较高。此外,文献显示研究区沉降速率一般为4~23.3 mm/a^[3-7]。我们利用PSInSAR技术提取到的沉降区位置分布和沉降速率跟文献显示的结果相吻合。

显示原图| 下载原图ZIP| 生成PPT

图6   PSInSAR与水准监测结果比较

Fig. 6   The monitoring results comparison of PSInSAR and leveling points

5 结 论

本文利用PSInSAR技术对现代黄河三角洲地面沉降进行监测研究,得到东城区和西城区平均地面沉降速率分别为7.4 mm/a和7.1 mm/a,沿海地区地面沉降速率约为4.4 mm/a,最大沉降速率达到22.1 mm/a,监测结果与水准观测数据吻合。通过与居民区分布图和油田分布图叠加分析,证实由现代黄河三角洲新成陆和冲淤特性决定的,地面沉降速率与地表载荷增加、石油开采活动和沉积物的自然固结压实具有较强的相关性。试验证明基于小数据集的PSInSAR技术可对湿地广泛分布的现代黄河三角洲地面沉降进行准确监测,此方法对环境复杂地区的地面沉降监测具有普遍推广意义。

致 谢：感谢Hooper提供StaMPS软件,荷兰Delft理工大学提供的Doris 软件及精密轨道数据,感谢USGS提供的SRTM数据。

The authors have declared that no competing interests exist.

参考文献

文献选项

• 原文顺序
• 文献年度倒序
• 文中引用次数倒序
• 被引期刊影响因子