基于GNSS信号的随机子空间模态参数识别方法

2020-12-22廖聿宸宗周红

东南大学学报（自然科学版） 2020年6期

吴睿廖聿宸宗周红张坤

(1东南大学爆炸安全防护教育部工程研究中心，南京 211189)(2中交第二公路勘察设计研究院有限公司, 武汉430056)

全球卫星导航系统(GNSS)是一种有效的结构变形监测手段，凭借自动化程度高、监测范围广、实时性良好等优势在桥梁健康监测领域得到推广应用[1-6].伴随着我国北斗卫星系统的完善及其在桥梁健康监测中的应用发展，基于GNSS信号识别桥梁结构模态特性的需求日益迫切.

GNSS信号受接收机噪声、电离层折射、对流层折射、多路径效应等多种误差源的干扰，且结构振动通常小于各种误差源的综合影响，因此需要对GNSS观测数据进行滤波处理，才能得到准确的结构振动信息[2].现阶段，常采用EMD及小波变换等方法对GNSS信号进行处理[7-11].

GNSS接收机能直接测量结构振动三维绝对位移，识别精度在毫米级，但受限于采样频率，不易识别高频振动.基于加速度数据识别模态频率和振型的方法相对成熟，识别精度较高，加速度信号经过2次积分可获得结构动态位移，但在积分过程中会产生趋势项，导致对低频成分不敏感，难以监测0.2 Hz以下的振动，丧失部分位移精度[2].因此，在有效提升GNSS信号采样频率和去除高频噪声的基础上，应采用其他更为稳定及准确的时频域参数识别方法对GNSS信号进行处理，从而获取结构高精度位移和多阶模态特性.

鉴于此，本文提出了一种基于GNSS信号的随机子空间模态参数识别方法.该方法利用小波变换在低频段频率分辨率高的特点[12]，对GNSS信号进行分频段滤波去噪.采用数据驱动的随机子空间模态参数识别方法对处理后的信号进行模态参数计算，并借助稳定图提取模态参数.最后，通过灌河大桥的模态分析对本文方法进行验证.

1 信号分频段降噪处理

小波分析[13-14]是将信号在不同的频域分辨率上进行分解，得到信号的近似分量和各分辨率下的细节分量，细节分量的数值即为小波系数[15].根据二进小波分解的特点，若信号采样频率为fs，对信号进行N层小波分解，则第k(0

小波阈值降噪技术将小于阈值的小波系数置零，再对处理后的小波系数进行重构，从而达到滤波降噪的目的[16].考虑到GNSS信号的采样频率高于可测的结构模态频率，若直接对其进行小波分解，需要采用较大的分解层数，而低层分量对应的频率范围不在关注的频段范围内，从而失去分析意义.因此，先对GNSS信号进行隔点取值的重采样，重采样频率不应小于关注频段最大频率的2倍.

GNSS信号滤波降噪处理步骤如下：

①设定关注的频段范围，对GNSS信号进行第1次小波分解，分解层次可根据最高层近似系数频段范围进行估计.取最高层近似系数进行重构，得到包含长周期位移及多种噪声的信号低频分量，去除该低频分量可保留大部分结构振动信息，提高关注频段范围内结构模态频率的识别效果[17].

②对高频分量再次进行小波分解，为减少细节分量中的噪声成分，可增大该次小波分解的层数[18].为保证频率识别效果，分析的频段范围应包含于关注的频段范围内，且选择分解层数较高的相邻2层小波系数进行重构，得到分频段重构分量.

③为适应2次小波分解后信号能量降低，采用改进的固定阈值对小波系数进行处理，即

(1)

式中，T为改进的固定阈值；σk为根据第k层小波系数计算所得的噪声标准差；Nk为第k层小波系数长度.

硬阈值函数处理后的重构信号与原始信号的逼近程度较好[19]，故选取硬阈值函数对小波系数进行处理.

2 数据驱动的随机子空间模态参数识别方法

数据驱动的随机子空间模态参数识别方法[20]将输出数据组成2mi×j的Hankel矩阵，借助QR分解，计算投影矩阵Oi，其奇异值分解为

(2)

式中，Oi为Hankel矩阵后i行对前i行的投影；U1∈Rmi×n、V1∈Rj×n为酉矩阵；S1∈R为对角阵，由非零奇异值组成；m为观测点数；i、j分别为Hankel矩阵的块行数和块列数；n为状态空间模型的系统阶数.

式中，Ψ为特征向量；ω角频率；λ为特征值；Λ为特征值λ组成的对角阵；t为离散信号的采样间隔.

对投影矩阵奇异值分解时，选取不同的奇异值数目，可计算出不同系统阶数的模态参数.稳定图方法[21]将不同系统阶数的模态识别结果绘制于同一幅图中，认为出现频次较高的稳定模态具有较大的可信度，并根据结构频率f、阻尼比ζ和实模态向量φ三类模态参数来评价模态的稳定性，假设状态空间模型的系统阶数范围为[nmin,nmax].

频率f采用分辨率判别，以0.01 Hz为区间长度将频段范围分成多个频率列，构造频率矩阵F.将各阶次频率计算结果存入相应的频率列中，当某一频率列内频率数目大于αNc时，认为该列频率稳定.其中，Nc为计算次数；α∈(0,1)为系数，由系统特点确定.

阻尼比存在空间和时间变异性，可用于剔除部分虚假模态[22].构造阻尼矩阵Z，将各阶次计算所得的阻尼比与阻尼阈值ζmax进行比较，将小于该阈值的阻尼比存入阻尼矩阵中.由于实际测量中存在测试噪声，且阻尼比对噪声比较敏感，在频率变化不大的情况下，阻尼比可能发生较大改变[22].因此，需对阻尼矩阵进行如下处理：

(3)

进行频率及阻尼比判别后，剔除虚假模态[23].构造二维矩阵J，其中Jpq=fpqζpq，Jpq、fpq、ζpq分别为矩阵J、F、Z中第p行q列的元素.假定结构满足Rayleigh阻尼形式，即

C1=aM+bK

(4)

式中，M为质量矩阵；K为刚度矩阵；比例系数a、b可用任意已知的二阶振型阻尼比确定，且b>0.

则有

(5)

式中，ωpq=2πfpq.因此，当系统阶数为p时，结构某二阶固有频率fpq、fpr应满足

fpq

(7)

将F、Z矩阵中对应元素置零，从而剔除各阶次部分虚假模态，使识别结果更加合理.

3 基于GNSS信号的灌河大桥模态参数识别

灌河大桥位于江苏省连云港市和盐城市交界处，主桥采用双塔双索面半漂浮5跨连续组合梁斜拉桥，主跨长340 m，立面图见图1(a).其健康监测系统于2013年底建设完成并投入运营，该系统测点包括三向加速度计、GPS接收机等传感器.其中，GPS接收机共设4个测站(见图1(b)中G1～G4)，包括主塔塔顶移动站2个、主梁跨中左右两侧各1个，采样频率为12 Hz，设参考站1个，布置在附近收费站处.三向加速度传感器共设11个(见图1(c)中A1～A11，)，采样频率为50 Hz.

(a) 立面图

(b) GPS测点布置图

选取测点G3、G4的竖向位移数据，数据时段为2020-02-18T00：30—01：00.采用本文方法对该段GNSS信号进行滤波处理和模态参数识别.

首先，设定关注的频段范围为0.30～1.10 Hz.选取sym3小波基函数对原信号进行3层小波分解，取第3层近似系数进行重构，得到低频分量，其频段范围为0～0.75 Hz.将原始信号减去低频分量后可得信号高频分量(见图2).由图可知，低频分量较原始信号更为平滑.对信号低频及高频分量进行频谱分析可知，分离低频分量后，余下的高频分量包含了结构振动响应及随机噪声，需对其进一步分解(见图3).

(a) 原始信号

(b) 低频分量

(a) 低频分量频谱

(b) 高频分量频谱

其次，选取sym7小波基函数对高频分量进行第2次小波分解，小波分解层数为5层.分别取第5层和第4层(D5+D4)、第4层和第3层(D4+D3)小波系数进行重构，信号频段范围分别为0.187 5～0.75 Hz、0.375～1.5 Hz.以第5层和第4层小波系数进行重构为例，将重构信号重采样为2 Hz的信号，结果见图4.

然后，选取sym12小波基函数对图5所示的信号进行3层小波分解，经小波阈值降噪后重构信号，结果见图5.

采用功率谱分析和数据驱动的随机子空间模态参数识别方法对处理后信号进行结构模态参数识别.图6给出了功率谱曲线和稳定图.取α=0.2，ζmax=0.1[24]，β=0.4[20]，若剔除虚假模态后的完全稳定点个数超过频率稳定点个数的30%，则输出该模态频率.

(a) 重采样信号

(b) 重采样信号频谱

(a) D5+D4阈值处理信号

(b) D5+D4阈值处理信号频谱

(a) D5+D4

(b) D4+D3

GNSS信号经分频段小波阈值降噪后，频段内噪声得到有效抑制，功率谱峰值更为显著.借助稳定图可识别出灌河大桥多阶模态频率，不同频段信号识别结果可相互验证；稳定图与功率谱曲线相结合，可提高识别结果的准确性.针对在某些结构频率附近存在多列稳定点的现象，比较各列稳定点所对应的阻尼比均值，若结果较为接近，认为上述稳定点表示同一模态，选择上述稳定点的频率均值作为该阶模态的频率.处理后的GNSS数据模态频率分析结果见表1.

表1 GNSS数据模态频率 Hz

4 基于加速度信号的灌河大桥模态参数识别

选取同一时间段、同一断面的竖向加速度数据进行模态参数识别.对加速度数据进行重采样，采样频率降至2.78 Hz，并进行频谱分析，结果见图7.由图可知，加速度信号噪声含量较少，在关注的频段范围内可观察到多个清晰峰值.

(a) 加速度重采样信号

(b) 加速度重采样信号频谱

图8给出了加速度信号的功率谱曲线与稳定图，各参数设置与第3节相同.灌河大桥的频率识别结果列于表2.其中，有限元模型的模态分析结果经过了环境振动试验修正[25].对比GNSS信号及加速度信号识别结果可以得出，借助2类信号均能识别出灌河大桥的模态频率，最大误差在3.44%以内，说明基于GNSS信号的随机子空间模态参数识别方法能够准确捕获斜拉桥的模态参数.