分类：农业论文 时间：2021-09-04 热度：202
提要:近年来，永久散射体合成孔径雷达干涉测量技术(PS-InSAR)快速发展，为宏观尺度的地表形变监测提供了契机，但该技术在国内尚未广泛应用，尤其在山体动态监测领域仍处于起步和探索阶段。文中利用PS-InSAR方法监测内蒙古和林格尔县的山体隆升情况，并对方法进行面向对象的改进。基于主被动遥感数据协同增加信息维度，构建多尺度滤波方法剥离人类活动及气候变化引发的地表形变，使监测对象锚定于山体隆升形变;再使用分水岭方法使监测窗口缩小至山脊区域，进一步提高了信噪比。监测结果表明，本区现今地质构造运动正处于差异性隆升阶段，东南部蛮汗山区以正向隆升形变为主，西北部土默川平原以负向沉降形变为主。其中，指示山体隆升的形变点多达745812个，占累计形变点总数的73.99%，平均隆升速率为17.59mm/a，构造运动特征为连续性应变、蠕变。该区近年来构造运动异常活跃，应在地震、地质灾害预测预警、重要工程选址时予以重视。
关键词:永久散射体时序差分干涉测量;山体隆升;多尺度滤波;动态监测;灾害预警
监测山体隆升变化对于研究地壳活动强度[1]、应力作用大小[2]、现今构造活动演化趋势[3-5]、地质历史溯源[6]、地震、地质灾害预测预警[7]、重要工程选址[8]等具有直观指示意义。目前，山体地质历史时期抬升速度研究主要依据综合古地理研究法和同位素年代学方法[9]。综合古地理研究法利用将今论古原则，通过对地层古生物、古土壤、地貌等多学科综合研究，恢复古地理环境，推测山体隆升幅度和时代。同位素年代学方法通过特征矿物或全岩的Rb-Sr、K-Ar或裂变径迹年龄，推算抬升速率及隆起年代[10]。两种方法均属于依据形变线索推测形变速度的间接监测方法。直接监测法主要基于精密水准测量[11]，需进行大量的地面测量工作，人力成本及经济成本较高，难以适用于复杂地形条件及大区域尺度，且测点密度稀疏，复测周期长[12]，难以满足精准监测山体隆升变化的组网密度和观测频数。
永久散射体合成孔径雷达干涉测量技术PS-InSAR(PersistentScattererSyntheticApertureRadarInterferometry)是一种精度更高、现势性更强的主动源遥感监测技术[13]，目前广泛用于资源调查[14]、环境监测[15]、灾害预报[16]和地形测绘[17]等领域。该技术通过分析时间序列影像，提取时序内相位和幅度保持稳定的离散点，即PS(PersistentScatterer，永久散射体)点，反演地表缓慢形变。通过长时间序列点目标分析技术，有效克服了空间去相干、时间去相干、大气延迟等系列问题，显著提高了测量精度和信噪比，高程向测量精度优于厘米级甚至达到毫米级[18]，对于监测宏观尺度的山体隆升变化具有较大优势。
当前，采用PS-InSAR技术监测山体隆升的研究案例较少，且面临诸如易混形变效应区分[5]、特征标识库建立、抬升机理解析等一系列理论及技术困难。文中拟利用PS-InSAR技术对内蒙古和林格尔县域内山体进行宏观动态监测，通过主被动遥感数据协同增加信息维度，构建时间域和空间域的多尺度滤波方法分离易混形变效应，并结合地质背景与地震数据对监测结果进行初步解释。
1材料与方法
1.1数据概况
采用主被动遥感数据协同方式监测内蒙古和林格尔县的山体隆升运动。主动遥感数据为2017至2018年的17期Sentinel-1A雷达数据，用于差分干涉以监测地体累计形变。雷达遥感数据选择干涉宽模式(InterferometricWide，IW)下的单视SLC(singlelookcomplex)数据，为C波段的升轨数据，C波段属于长波，可以穿透稀疏植被直接作用于地表。该数据处理后的距离向分辨率为2.3m，方位向分辨率为13.9m，极化模式采用VV极化。被动遥感数据为同步观测的Sentinel-2A数据(10m)。辅助数据为DEM数据(30m)、土地利用数据(1:1万)、地质背景数据及地震数据(50a)。
1.2研究方法
1.2.1技术路线
PS-InSAR方法监测山体现今隆升变化的主要技术流程是:1)主被动遥感数据协同处理，包括数据导入、精轨校正、主影像筛选、几何配准、子条带拼接，光学遥感数据反射率生产、拼接镶嵌、联合配准等。2)差分干涉，包括设置时空基线阈值、求解DEM误差、迭代更新初始DEM、去平地相位和去地形相位、生成差分干涉图等。3)反演永久散射体累计形变，包括提取永久散射体(即PS点)、构建Delaunay三角网、计算PS点的参数信息、求解二次差分相位模型和PS点绝对形变量、残余相位解缠、大气滤波、线性形变量和非线性形变量计算、地理编码等[19]。4)多尺度滤波，包括空间域滤波及时间域滤波，排除人类扰动和外部自然变化引起的地表形变，将ROI锚定在山体隆升形变。5)综合解译，对多源数据进行信息复合和KDG(knowledgediscoveryfromGIS)分析[20]。总体技术路线见图1。
1.2.2差分干涉处理
主从影像共生成17期差分干涉图，依次执行感兴趣区选择、干涉处理、地形相位去除、差分模型构建四个步骤。主辅影像干涉相位可表示为公式(1)。式(1)中φtopo为地形相位;φdef为形变相位;φatm为大气延迟相位;φflat为平地相位;φnoise为噪声相位。在差分干涉基础上，利用Delauany三角化方法构建三角剖分网络，基于PS点的三角剖分网，计算网络上每条边的相邻PS点的二次差分相位;然后进行相位解缠，得到每个PS点的绝对形变速率和DEM修正量，进而反演累计形变量;最后使用地理编码将形变解算结果从雷达极坐标系统转换到地理坐标系统。φ=φtopo+φdef+φatm+φflat+φnoise(1)
1.2.3多尺度滤波
地表累计形变可能由多种因素构成，包括人类扰动和自然变化，自然变化又包含土壤侵蚀、滑坡、泥石流、冻土胀缩、岩石风化、搬运、堆积等。因此，需要对环境噪声进行多尺度滤波，剥离人类扰动和自然变化影响，将监测对象锚定为因地球内部活动引发的地表形变。
多尺度滤波的滤波方法主要为空间域滤波和时间域滤波。空间域滤波又分为植被滤波、稳定地表滤波、人类活动滤波、堆积区滤波。首先，通过计算主末影像的相位标准差，加入植被指数的限制因素，振幅小的信号通过，振幅大的剔除，滤除茂密植被，将感兴趣区限制于植被稀疏区域。其次，设置稳定地表值域，对永久散射体的累计形变信号进行高通滤波，滤除无变化或变化极小的区域。然后，根据土地利用数据掩膜城市村镇的不透水面及活跃耕地，并设置阈值分割上限，对永久散射体的形变信号进行低通滤波，滤除人类活动干扰。最后，使用分水岭方法进行堆积区滤波，依次进行山谷填洼(Fill)、流向计算(Flowdirection)、流量计算(Accumulation)、累积流量分割(Gridcalculator)、分水岭(Connection)等处理，分离出山谷洼地等沉积堆积区域。
2结果与分析
2.1山体隆升形变点筛选结果
对2017至2018年的17期雷达遥感影像进行时序差分干涉，共计提取出1007958个累计形变点，覆盖县域全境。每个(或数个)累计形变点对应一处形变区域。地表形变场为多重驱动力下的复杂形变场，从数以百万计的形变点中筛选出反映山体隆升变化的形变点，需要一个再滤波的过程。对所有累计形变点执行时空域多尺度滤波，识别出反映山体隆升变化的形变点。空间域滤波依次执行植被滤波、稳定地表滤波、人类活动滤波、堆积区滤波。植被滤波基于相位标准差低通滤波滤除高植被覆盖区域，主要为活跃耕地及森林。稳定地表滤波将年累计形变量小于±1cm区域视作稳定地表[21]，将变化量小于此数值的累计形变点剔除。人类活动滤波基于土地利用分类结果，将全部人类活动区域剔除，如建筑用地、设施用地、公路铁路、耕地、镇村等，只保留自然地表。堆积区滤波基于DEM进行山谷填洼处理，避免与山谷区域的沉积堆积效应引发混淆。空间域滤波结果见图2。
时间域非线性滤波基于累计形变点在时间序列上的变化特征，分离出季节性、气候性形变点。畸变滤波需设置畸变阈值，随机抽取正向运动且年累计形变量大于1cm/a的永久散射体样本5612个，分析样本库数值统计特征。样本平均形变量在40mm/a后迅速收敛接近零值，故设置畸变阈值为40mm/a，超过此阈值的形变点作为畸变点剔除。经由多尺度滤波处理，将大部分人类扰动或自然变化引起的地表形变点视作噪点剔除。经筛选，共提取指示地体隆升形变的累计形变点合计745812个，占全部累计形变点总数的73.99%，平均隆升速率为17.59mm/a。
2.2山体隆升形变时空特征
山体隆升形变点与高分光学遥感影像及坡向图叠加分析可以划分出2个主要构造地貌单元(图3):县域西北部的土默川冲积平原区，约占县域面积的22%，由黄河、黑河冲积而成，地形平缓;县域东南部的黄土丘陵-蛮汗山区，约占县域面积的78%，相对高差200-300m，绝大部分山体坡度小于35°。
本区地体抬升的空间差异明显，累积抬升点云主要分布在黄土丘陵-蛮汗山区，极少分布在土默川平原区。一方面平原大部分地区的人类活动密集，人类活动区域(不透水面、活跃耕地等)已由多尺度滤波环节去除;另一方面少数具有观测条件的平原区的累计形变多为沉降变化，为构造运动、土壤侵蚀和地下水超采等因素综合所致。监测结果显示，本区东南部的黄土丘陵-蛮汗山区主要为隆起抬升的正向运动，西北部的前陆平原主要为沉降凹陷的负向运动。绝大多数隆升形变点沿山脊线呈叶脉状分布(图3)，一方面可能与局部构造应力方向相关，山脊形态主要受应力方向控制，高程形变显著;另一方面山脊区域基岩广泛出露，相干条件较好，识别出的永久散射体数量较多。
相关期刊推荐：《干旱区资源与环境》(月刊)1987年创刊，是综合性学术刊物。凡研究干旱半干旱甚至季节性干旱区的论文，研究干旱半干旱问题及防治技术的论文，特别是绿洲建设的论文，不论属于社会科学领域或者自然科学领域均可刊用
雷达时序差分干涉通常需要大于15期影像以建立联合方程并迭代求解，为捕捉地表形变的季节规律和物候特征，也需要时间序列上的高采样频率。本研究参与差分干涉的影像为17期，时间分辨率优于月度。对抬升形变点云的月度变化量取平均值，回溯抬升点云在时间序列上的变化。由图3柱状图可知，在12个月中，山体累积抬升量逐月均匀增加，不受季节性气候变化及人类扰动影响，反映为地球内部构造应力驱动，总体呈线性增加趋势。根据均衡理论[22]，地表形变的主要驱动因素是构造隆升和侵蚀夷平，山体构造隆升后，受剥蚀减荷作用而自重减少，在基底应力推动下会继续隆升，因而构造抬升通常与均衡抬升并存，在时间上表现为连续性、渐进性形变。总体而言，研究区境内山体正处于连续性、差异性隆升阶段。
3讨论
3.1监测结果可靠性论证
雷达遥感适用于监测高程形变信息，光学遥感适用于监测平面纹理信息，二者结合可以实现景观多要素三维动态监测。山体隆升监测需要将景观多要素形变信号分离，融合使用多源数据并提高监测频度是合理的选择。针对山体现今隆升运动的一般性时空特征，文中对方法应用进行了一些改进，通过多尺度滤波剥离了人类扰动和外部自然变化导致的地表形变效应，使监测对象锚定于山体隆升形变，计算得出本区现今山体隆升速率为17.59mm/a。
InSAR方法精度验证主要是形变径迹调查[23]，和与高精度GPS或精密水准测量数据进行比对[4]，课题组于2018年5月、2019年8月两次赴研究区实地调研，采用高精度GPS测量、形变径迹追索、实地拍照等方式进行实地采样103点，GPS高程测量精确到2位小数(厘米级)。经实地验证，86.4%的形变点可以追索到诸如滑坡、侵蚀等形变径迹，总体测量精度控制在厘米级，说明所述方法的监测结果具有可靠性。
3.2山体隆升成因浅析
研究区大地构造位于华北板块北缘，燕山-阴山板内造山带西段，其特征主要与燕山-阴山造山带有关。断裂是研究区的主要构造，块断结构是基本的构造格架。NE、NNE、NW向断裂控制了整个研究区的构造格局，将研究区分为后公喇嘛中新生代断陷盆地、和林格尔中新生代断陷盆地、胜利营-永兴乡古隆起三个构造单元[24-25]。其中，NE向黑山泉-西门沟高角度正断层将研究区分为两个不同的地质地貌单元，北西侧的平原区和南东侧为丘陵-山区，与山体隆升形变点与高分光学遥感影像及坡向图叠加分析结果一致。结合地质背景资料分析，研究区的山体隆升可以追溯到燕山运动时期[26]，晚白垩世末期燕山运动使全区隆起，进入新生代早期阶段区内一直遭受剥蚀和夷平，形成了新近纪古夷平面。中新世火山活动异常强烈，研究区喷发了大面积的汉诺坝玄武岩。进入第四纪以来，山区以间歇性抬升为主，形成了多级夷平面和黄土丘陵等构造地貌单元。晚更新世之后地壳运动渐缓，黄土丘陵区堆积了马兰期风成黄土，此时古气候环境已由温暖潮湿转变为干燥强季风环境。全新世以河流相沉积为主，但差异性升降依然具周期性变化，全新世以来山体不断隆升，从而形成了山前广为发育的冲洪积平原和不同的河流阶地与河床沉积。综合文中监测结果，本区地体形变在空间上呈差异性特征，东南山区隆起，西北平原沉降;时间上呈连续性渐变特征，年均隆升速率较快，直接显示该区近年来构造活动异常活跃。此外，中国地震台网地震监测数据显示，研究区自20世纪70年代以来共发生46次地震(1次灾害性地震，1976年4月6日和林格尔6.2级)，其中2018年发生中-小震24次，进一步说明本区现今构造运动正处于异常活跃期，断裂交汇部位、转折端是地震易发点，需高度重视[27]。
4结论
文中以星载合成孔径雷达差分干涉测量为主要技术手段，通过主被动遥感协同监测、多源数据信息复合、多尺度滤波等方式，构建了在宏观尺度上实现山体隆升快速监测的技术框架。以内蒙古自治区和林格尔县为案例区开展实证研究，通过实地调查对方法可靠性进行验证，并结合地质背景资料，对监测结果进行了初步解释。研究得出结论如下:
(1)内蒙古和林格尔县2017-2018年山体平均隆升速率为17.59mm/a，主要为构造运动引起，地体隆升幅度较高，是华北北部中低山地差异性隆升的典型地区，具有代表性意义。
(2)和林格尔县区域构造运动方式以连续性应变、蠕变为主，小震频繁，现今地质构造运动异常活跃，建议在重大项目建设选址、地震及地质灾害预测预警中予以重视。——论文作者：杨通1，南燕云2，张智杰3，于潇4，王效科1
文章名称：基于主被动遥感协同的山体隆升动态监测研究
