基于低成本毫米级GNSS技术的膨胀土边坡现场监测

2022-02-23陈孜黄观文白正伟张双成张永志

中南大学学报（自然科学版） 2022年1期

陈孜，黄观文，白正伟，张双成，张永志

(长安大学地质工程与测绘学院，陕西西安，710054)

我国膨胀土面积约占陆地总面积的1/3，每年因膨胀土工程地质灾害造成的经济损失达数百亿元[1]。膨胀土是一种富含强亲水性的蒙脱石黏土矿物的土体，易在干湿循环作用下表现出吸水软化、失水收缩开裂的反复涨缩形变特征[2-4]。膨胀土变形特点决定了膨胀土滑坡造成的破坏具有反复性、碎片化性和长期潜在性[5]。由膨胀土滑坡的浅层性特点可知，地表位移是监测的关键，对具有隐患的膨胀土滑坡开展实时地表位移监测是掌握其位移运动和演化过程的前提。

国内外已有许多学者对膨胀土滑坡监测进行了相关研究。缪伟等[6]利用测斜仪对某膨胀土边坡进行了现场变形观测，研究了膨胀土边坡变形随季节演化的过程。张新生[7]利用测斜仪和PVC测斜管对膨胀土深部位移进行监测分析，确定了滑动面位置。这些常用的监测方法虽然在确定滑动面及深部位移信息中有独特优势，但也存在不足。首先，测斜仪或测斜管不能监测坡体的三维形变信息。其次，施工效率较低且监测成本较高[8]。对于安全隐患区面积较大的膨胀土滑坡，难以实现大范围监测点布设。

全球卫星导航系统(GNSS)技术可提供实时全天候三维定位，已经逐渐成为滑坡监测的常用方法。随着我国北斗三代系统迅速完成全球组网，自带短报文通讯功能的北斗卫星导航系统(BDS)在滑坡监测预警中发挥了越来越重要的作用。相比全站仪等其他滑坡监测技术，GNSS尤其是BDS具有自动化、实时、高精度等显著优势[9-12]。但GNSS 应用于膨胀土监测也存在一定的技术局限性。首先，传统GNSS 实时精密单点定位(PPP)和实时动态差分定位(RTK)技术精度为厘米级[13-14]，在毫米级要求的膨胀土监测领域应用尚有困难。其次，膨胀土滑坡分布不集中，需要大范围布设监测设备[15]，因此要求监测设备体积小、成本低；而普通GNSS监测设备以立杆式监测形态为主，空间体积较大且成本昂贵。第三，膨胀土滑坡周边易受到植被和建筑物遮挡，多路径误差及环境误差干扰大，影响GNSS定位精度。

针对上述难题，本文作者设计一种适用于膨胀土滑坡的低成本GNSS监测装备，并在监测立柱中加入可变换高度的天线，设计多测站共享电源等，以实现监测站低成本和小体积的目标。针对环境对信号的干扰问题，提出一种适用于膨胀土监测的抗干扰数据处理策略，通过实时对监测误差进行修正，保障精度。最后，以安康机场膨胀土试验区为例，验证GNSS技术在膨胀土滑坡长期实时监测中的适用性和精度。

1 GNSS监测的抗干扰技术

1.1 GNSS监测终端及装置

针对传统GNSS监测设备体积大、成本高的问题，本文提出“简化监测设备+共享电源”的思路，从设备体积、成本和功耗3个方面入手，研发出一种适用于膨胀土地表位移监测的低成本GNSS监测设备，该设备由小型化GNSS监测装置和共享电源集成站2个部分组成。

小型化GNSS 监测装置由GNSS 监测终端(接收机)、信号接收装置(天线)、连接装置和侧边反射棱镜构成。设备结构设计见图1。装置的核心是GNSS 监测终端(接收机)，固定于监测装置底部。GNSS 监测终端中剥离了非必要功能模块(如定位解算、数据存储、电池等功能模块)，仅保留数据采集模块和通讯模块，通过数据采集单元收集各监测点的实时观测数据，通过GPRS/4G/5G等无线通信手段将监测数据实时传输至北斗云监测平台，基于“传感器+云端服务”的理念，利用北斗云在线解算取代监测终端的自解算模块，实现云存储和云解算。此设计不仅使得终端体积缩小了75%，同时也有效地降低了硬件成本。自研设备与不同品牌接收机的参数对比见表1。除此之外，为了保证成功捕获滑坡初始阶段到加速阶段的完整形变数据及降低监测功耗，GNSS监测终端加入了自适应变频功能。在膨胀土滑坡运动处于稳定状态时，以较低的监测频率(如1 次/h)进行常态化监测；而当坡体遭遇降雨或气候变化，膨胀土体易处于非稳定状态时[15]，GNSS监测终端自适应地加快采样频率，例如2 次/s。自适应变频技术不仅可以捕获更加全面的监测数据，还可以降低长期不间断监测的功耗。

图1 监测装置设计Fig.1 Design of monitoring device

除了缩小BDS 监测终端的尺寸外，本文监测装置将传统的立杆式监测杆(高度1.8 m 以上)替换为可调整高度(10～50 cm)的小型天线连接杆。根据不同边坡观测环境，自适应地调节天线连接杆高度，这在很大程度上减少了公共环境干扰问题。为了无缝匹配测量机器人进行同期监测数据检核，在天线连接杆处侧边安装测量机器人棱镜贴片，棱镜贴片可水平旋转360°，上下旋转30°。

在设备供电方面，采用立体式“共享电源”集成站设计，见图2。共享电源集成站的主要作用在于向膨胀土监测区内各GNSS监测终端提供长期稳定供电。设计思路如下：将太阳能蓄电池与GNSS监测终端/GNSS天线装置剥离，并在坡体后缘安放“一对多”集成电源站，通过电缆线与边坡监测装置连接。一个共享电源站可以对2～3个小型GNSS监测装置长时间稳定供电。这不仅减少了监测装置的体积，而且是集柔性太阳能电池技术、高性能胶体蓄电池于一体的高效率供电解决方案。在滑坡体后缘稳定处建立共享电源集成站，单个共享电源集成站可同时向多个GNSS监测点提供连续稳定供电，并根据实测区域调整阳光照射的角度以获取足够多的电力资源，保证了该监测设备的长期低功耗监测和稳定供电。

图2 共享电源集成站组成Fig.2 Composition of shared power integration station

1.2 高精度实时监测抗干扰技术

将上述低成本GNSS监测设备应用于膨胀土滑坡监测中会遇到3 个技术问题。首先是精度问题，传统RTK/PPP 技术只能实现厘米级定位精度，且遇到通讯不稳时监测结果常常不连续。膨胀土滑坡短期变形小，一般仅为毫米级且须连续稳定监测，因此必须提升监测定位算法的精度和稳定性。其次是环境遮挡问题。膨胀土滑坡周边常有茂密植被或建筑物遮挡，会造成监测定位结果中存在毫米至厘米级的系统误差。第三是信号反射问题。监测设备的小型化导致接收天线距离地面较近，会明显受到多路径误差的影响[16-17]。

针对上述3个技术问题，本文基于自研小型化GNSS膨胀土滑坡监测设备，提出适用于膨胀土监测的GNSS 抗干扰数据处理策略。主要设计思路为：1)监测数据实时传输及在线质量分析；2)监测站多路径误差建模；3)测站环境建模补偿；4)GNSS实时监测质量控制。

1.2.1 监测数据实时传输及在线质量分析

膨胀土GNSS监测设备实时接收测区观测数据与卫星星历数据，基于Ntrip 传输协议将实时数据流传输至北斗云平台。北斗云平台对实时数据进行解码，从数据完整率、位置精度因子SPDOP、多路径效应等多指标对GNSS观测数据进行在线数据质量分析，剔除精度较差的伪距和载波观测值，并屏蔽不可用卫星等。

多路径效应主要反映监测站周围反射物体对接收GNSS信号的干扰误差。为了评估监测站受到的多路径误差，需进行多路径效应的计算。在实时定位应用中，通常忽略载波相位多路径误差，只考虑伪距多路径误差的影响，其计算公式为：

式中：fi和fj为频点i和j对应的频率；Mi和Mj为频点i和j对应的伪距多路径误差；pi和pj分别为频点i和j的伪距观测值，φi和φj分别为频点i和j的载波观测值。在监测区采用双频观测值，数据处理中一般将M1的阈值设为50 cm，将M2的阈值设为75 cm。

数据完整率R反映的是信号实际观测数量F与理论观测数量N之比，其计算公式为

数据处理中，R的阈值设置为80%。数据完整率除了可以评估接收机性能之外，也可以表征监测站周边的环境优良程度。将超过阈值的观测值进行剔除处理。

GNSS位置精度因子SPDOP可衡量监测站上空卫星星座分布状况，其计算公式为

式中：σXu，σYu和σZu分别为测站三维坐标的均方根误差分量；σuere为伪距误差。以上评估指标均可以分析单测站周边环境及天空星座分布状况。

1.2.2 监测站多路径误差建模

不同监测站间的多路径误差没有相关性，GNSS差分定位不能消除多路径误差，小型化监测设备易受到多路径影响。本文提出基于多天观测数据，利用多路径误差周期重复性特点[18-21]，采用经验模态分解(empirical model decomposition，EMD)或小波变换等方法提取前一个回归周期多路径误差，并将提取结果应用到下一周期以实时改正动态定位，减少多路径误差对膨胀土实时监测的影响。

经验模态分解(EMD)算法具有良好的时频分辨率和适应性，能够很好地重建原始信号。EMD 算法可以将输入信号分解成相应的高频和低频分量，而多路径误差主要表现为低频特征。通过选取适当阶数对信号进行分解并重构就可以获得多路径误差及噪声影响。信号x(t)分解模式为

其中：n为分解尺度；为残差序列；为第s个模态函数；t为时间序列中的某一时刻。

信号分解步骤如下：

1)找出x(t)的所有极值点；

2)用插值法形成上下2 个包络函数emin(t)和emax(t)；

3)计算均值m(t)=[emin(t)+emax(t)]/2；

4)抽离细节信息d(t)=x(t)-m(t)；

5)对剩余的均值重复步骤1)至步骤4)。

1.2.3 测站环境建模补偿

膨胀土滑坡监测环境易受植被、山体等环境遮挡等干扰，若依据传统预置高度截止角将导致监测解算精度降低，此时需要构建基于周边地形的自适应高度截止角Em。可利用监测站与可见卫星的高度角、方位角等信息对监测站周围空间地形环境进行自适应建模[22]，将建模信息存入北斗云平台。采用适合多遮挡情况下的部分模糊度固定策略提高计算效率。

基于测站周边地形空间所建模型的表达式如下：

式中：Em为基于测站周边空间地形的自适应卫星高度截止角；Ec为默认高度截止角，通常为常量；Ea为表征方位角变化的地形特征信息。通过对足够多的观测离散量进行曲线拟合即可得到Ea。自适应高度截止角由周边地形线与监测站形成的高度截止角加上一个浮动角度构成，可提高复杂监测环境下的解算精度。

1.2.4 GNSS实时监测质量控制

根据Baarda 数据探测理论，将信号观测值进行相应线性组合，对原始数据进行粗差探测和标定[23]，即可实现对实时观测数据的质量控制。根据不同基线长度和测站间高差估计电离层误差与对流层误差，优化函数模型。根据环境建模信息自适应优化随机模型。利用PPP技术对基准站进行稳定性分析和基准补偿，并利用基线长度及大气约束信息进行GNSS动态增强定位。滤波后进行验后残差检验，提高监测序列的可靠性。

GNSS滑坡监测质量控制主要包括前端质量控制与后端质量控制，其中多路径的抑制属于前端质量控制，后端质量控制主要包括信号观测值粗差的检验和剔除。粗差探测基于均值漂移模型，该模型认为含粗差的信号观测值的期望发生变化，方差不变，其残差服从正态分布。假设第l个信号观测值含粗差，则相应的均值漂移模型为

式中：L为信号观测值；r为漂移量；ε为观测噪声与未模型化噪声误差，其方差为σ2；B为列满秩设计矩阵，X为参数向量。利用Baarda数据探测法对信号观测值进行线性拟合，其原假设H0与备用假设H1如下：

在原假设条件下构建标准化残差=(其中，vl为信号观测值残差，σ0为信号观测值的标准差，Cll为残差的方差-协方差矩阵的第l个对角线元素)，当＞α(显著性水平)时，则认为信号观测值含有粗差，拒绝原假设H0，据此剔除信号观测值中的粗差。保留模糊度固定后的解算结果并对其进行卡尔曼滤波平滑处理以提高解算精度。

2 安康膨胀土边坡监测分析

2.1 工程概况

安康膨胀土边坡位于陕西省安康市汉滨区丁家沟机场快速干道南侧。该边坡位置及全貌分别如图3和图4所示，该区域附近膨胀土分布较广，边坡上下高差约33.4 m，边坡倾角约46°。经现场地质勘察可知该边坡土体为中膨胀土。由于边坡坡度较陡，易对机场干道造成滑坡隐患，已对该边坡采取了相关加固措施，但依然有变形滑坡的危险。本文作者对该膨胀土边坡布设1 个GNSS 基准站及8 个GNSS监测站进行实时监测，基准站与监测站分布分别图5和图6所示。监测点位置分布见图7。该边坡为经过削坡之后的治理边坡，分为上、中、下3层。边坡顶层布设2个监测点，边坡中层布设3个监测点，边坡下层同样布设3个监测点，监测点均位于各坡层前缘。

图3 安康膨胀土边坡位置Fig.3 Location of Ankang expansive soil slope

图4 安康膨胀土边坡全貌Fig.4 Overall view of Ankang expansive soil slope

图6 GNSS监测站Fig.6 GNSS monitoring station

图7 膨胀土边坡监测点位置分布Fig.7 Distribution of monitoring points on expansive soil slope

2.2 监测结果分析

2.2.1 实时监测精度分析

为了揭示膨胀土边坡GNSS实时监测精度，在不做平滑处理的前提下，计算监测点(P1)在单天内(2021-11-01)实时秒级采样下的GNSS三维变形量序列。以测量机器人监测结果作为真值，得到P1监测点1 d内的实时变形量残差即定位精度，如图8所示。统计3 个方向变形量残差的中误差eRMS即监测精度，结果如表2所示。

图8 GNSS三维方向实时定位精度Fig.8 Real time positioning accuracy of GNSS in 3D direction

表2 膨胀土滑坡GNSS监测精度Table 2 GNSS monitoring accuracy of expansive soil landslide mm

从图8可以看出，不平滑的原始GNSS实时秒级监测结果仍能够稳定保持在毫米级量级。由表2可见：GNSS在东方向、北方向、高程方向的单天实时定位精度分别为2.0，3.0和4.5 mm。由上述监测结果分析可以得知，当膨胀土边坡发生毫米级及更大形变时，本文BDS/GNSS 高精度监测技术能有效实时捕捉形变结果。

2.2.2 安康膨胀土边坡变形特征

GNSS实验数据采集时间为2020-11-01—2021-02-15。图9～10所示分别为膨胀土边坡监测点P1和P2在东、北及高程方向上的累计位移和位移速率。由图9(a)可知：采样时间段内P1监测点在东方向和北方向2 mm 内小幅变形，垂直方向位移变化幅度稍大，但稳定后其形变仍处于毫米级量级。由图9(b)可知：P1监测点整体位移速率较小且位移速率上下波动，在东、北方向变化幅值不超过0.5 mm/d，高程方向变化幅值不超过1.0 mm/d。P2监测点位移情况与P1的类似，位移速率曲线以0为中心波动且在东、北方向变化幅值不超过1.0 mm/d，高程方向变化幅值不超过2.0 mm/d。

图9 P1监测点累计位移与位移速率随时间的变化Fig.9 Cumulative displacement and displacement rate variations of P1 with time

图10 P2监测点累计位移与位移速率随时间的变化Fig.10 Cumulative displacement and displacement rate variations of P2 with time

黄土滑坡位移在东、北和高程3个方向上的变化具有一致上升或下降的特点，而膨胀土易受气候、水分和土质等影响，会产生周期性膨胀或收缩，这会导致监测点3 个方向的位移变化不一致[24]。本文P1和P2累计位移呈现出周期性震荡现象也反映了膨胀土的这一特殊变形特征，这与膨胀土体干湿循环条件的涨缩效应是相符的。

2.2.3 降雨对膨胀土变形的影响

膨胀土的涨缩变形与降雨关系密切。图11所示为安康P1与P2监测点垂向位移与降雨量的关系。从图11可以看出，每当该膨胀土边坡经历降雨之后，监测点的垂向位移都会出现上升趋势，且降雨对边坡形变的影响存在一定滞后性，这与膨胀土体本身遇水膨胀的规律吻合。当边坡处在非降雨时期时，垂向位移又表现出回落趋势，这与膨胀土失水收缩的特点相吻合。

图11 P1和P2监测点累计位移与降雨量的关系Fig.11 Relationship between cumulative displacement and rainfall of monitoring position P1 and P2

由于监测周期内处于旱季，降雨量较少。对3 d 累计降雨量与累计位移进行相关性分析，结果如图12所示。此次采用皮尔逊相关系数分析降雨与位移的滞后相关性。皮尔逊相关系数法很适合衡量2个标准差都不为0的连续变量(如降雨量与位移)的相关程度。

式中：ρU,S为变量U和S的皮尔逊相关系数；U和S为随机变量；Uˉ和Sˉ分别为U和S的期望。

从图12可以看出，与黄土滑坡类似，降雨对膨胀土位移变化存在滞后性，其中P1监测点垂向位移与降雨的相关性在滞后17 d 后达到峰值0.56，P2监测点垂向位移与降雨的相关性在滞后16 d后达到峰值0.57。对于安康机场膨胀土试验区，滞后周期为15～20 d，降雨量与形变为中度相关。该滞后周期比黄土滑坡的稍长，这与监测周期处于冬季，降雨量较少及干湿变化频率不高有关。

图12 P1和P2监测点位移与降雨量的滞后相关性Fig.12 Lag correlation between displacement and rainfall of monitoring position P1 and P2

2.2.4 相对湿度对膨胀土变形的影响

与降雨一样，相对湿度同样促进了土体的吸水膨胀效应，间接影响膨胀土变形趋势。对安康P1和P2这2个监测点相应时段的累计位移与相对湿度数据进行分析，结果见图13。由图13可见：在监测周期内，P1和P2变形趋势在10～70 d的监测时间段内有着明显的相同变化趋势，可以判断两者存在一定相关性。对两者结果进行滞后相关性分析，结果如图14所示。由图14可看出：当滞后周期为5～7 d时，实验区域内相对湿度对位移变化的影响达到峰值，P1的相关性极值为0.50，P2的相关性极值为0.53，为中度相关。