一种基于卡尔曼滤波数据融合的动态位移获取方法

摘要

本发明属于土木工程结构健康监测技术领域，并公开了一种基于卡尔曼滤波数据融合的动态位移获取方法。该方法包括：(a)将待处理对象简化为悬臂梁，采用多个节点将所述悬臂梁划分为多个子单元，构建悬臂梁的节点处位移与应变值与虚弯矩之间的关系式；(b)构建以某一节点处位移和速度作为变量的状态向量，建立状态向量预测值与上一时间步的状态向量修正值和加速度之间的关系式；(c)构建状态向量修正值关于状态向量预测值、修正系数和某一节点处位移的关系式；(d)实时采集悬臂梁的每个节点处对应的应变和加速度，计算获得状态向量修正值，以此获得最优动态位移。通过本发明，解决结构位移计算复杂和精度低的问题，抗噪性高，工程适用性强。

说明书

技术领域

本发明属于土木工程结构健康监测技术领域，更具体地，涉及一种基于卡尔曼滤波数据融合的动态位移获取方法。

背景技术

细长结构(超高层建筑、桥塔等)与普通结构相比，具有很高的高度，在强风和地震等荷载作用下会产生较大的水平位移，进而容易引起结构损坏或失稳，影响结构的可靠性和安全性。因此，对结构的水平位移监测与控制是超高层建筑、桥塔等大型结构健康监测的重要内容，精确的结构位移测量方法是进行位移监测和控制的关键。

位移测量的一般方法包括直接测量和间接测量方法，由于细长结构独特的结构体系、动力特性和周边环境，这些方法应用于大型细长结构时还存在较大局限性，主要体现在以下几方面：

(1)直接测量方法基于位移计、GPS、激光全站仪和机器视觉等设备，其中位移计需要安装固定支撑，GPS在大型结构的位移测量中精度不够高，激光全站仪和机器视觉方法容易受到天气因素和四周环境的影响，在大型结构位移测量中存在很大局限性。

(2)间接测量方法一般通过测量结构的加速度或应变来计算结构位移。加速度方法可计算得到结构动态高频位移，但容易产生基线偏移而使位移不准确；应变方法一般基于结构振型函数，获取大型结构的振型函数需要建立模型进行模态分析，比较耗时。

近年来，国内外学者尝试使用数据融合方法进行位移计算，也即将不同的位移测量方法结合，取长补短，从而提高位移结果的精度。Smyth等人提出了一种基于卡尔曼滤波的位移计算方法，将位移与加速度融合得到更精确的位移；Kim等人(2014)利用卡尔曼滤波将位移与速度融合计算高精度位移；对于细长结构而言，直接测量方法因存在上述局限性而不适用，因此考虑将不同的间接测量方法结合以计算结构位移。Park等人基于FIR滤波将加速度和应变结合计算简支梁桥的位移；Cho等人基于卡尔曼滤波将加速度和应变计算所得位移进行融合。上述数据融合方法均需要结构的振型函数，其精度受振型函数准确性影响。

发明内容

针对现有技术的以上缺陷或改进需求，本发明提供了一种基于卡尔曼滤波数据融合的动态位移获取方法，通过建立动态位移和加速度与动态位移之间的关系式，实现通过测量加速度和应变计算获得结构的动态位移，由此解决计算动态位移中需要考虑结构的振型函数，且精度受振型函数准确性影响的技术问题。

为实现上述目的，按照本发明，提供了一种基于卡尔曼滤波数据融合的动态位移获取方法，其特征在于，该方法包括下列步骤：

(a)将待处理对象简化为悬臂梁，采用多个节点将所述悬臂梁划分为多个子单元，构建第k时间步所述悬臂梁的节点处位移x_m(k)与应变值与虚弯矩之间的关系式(一)；

(b)构建以位移和速度作为状态变量的状态向量，利用所述状态向量将第k时间步的位移与第k-1时间步的位移和加速度之间的关系式(二)转化为第k时间步状态向量与第k-1时间步的状态向量和加速度之间的关系式(三)，根据该关系式(三)建立第k时间步所述状态向量预测值X(k)与第k-1时间步的状态向量和加速度之间的关系式(四)；

(c)根据所述状态向量计算该状态向量的卡尔曼增益，并以此作为所述状态向量预测值的修正系数K(k)，利用该修正系数构建第k时间步所述状态向量修正值关于所述状态向量预测值修正系数K(k)和节点处位移x_m(k)的关系式(五)；

(d)实时采集所述悬臂梁的每个节点第k时间步处对应的应变和加速度，利用所述关系式(一)～(五)计算获得第k时间步的状态向量修正值，从状态向量修正值中获得第k时间步的位移，以此实现动态位移的获取，其中，k是时间步，k＝1、2、3……。

进一步优选地，在步骤(a)中，所述关系式(一)优选采用下列表达式：

其中，i是节点的数量，p是第i个节点前节点的总数量，是第i个节点的虚弯矩，是第i-1个节点的虚弯矩，h_i是第i个子单元的长度，s_i是第个i子单元的悬臂梁的宽度，Δε_i(k)是第i个节点在第k时间步的应变，Δε_i-1(k)是第i-1个节点在第k时间步的应变。

进一步优选地，在步骤(b)中，所述关系式(二)优选按照下列表达式进行：

其中，是状态转移矩阵，是控制项输入矩阵，是系统噪声矩阵，Δt是预设的k与k-1时间步之间的时间间隔，w_d是位移系统误差，w_v是速度系统误差，X(k)是k时间步的状态向量，X(k-1)是k-1时间步的状态向量，是测量获得的第k-1时间步的加速度。

进一步优选地，在步骤(b)中，所述关系式(三)优选按照下列表达式进行：

其中，是k时间步的状态向量的预测值，是k-1时间步的状态向量的修正值，是在k-1时间步测量获得的加速度。

进一步优选地，在步骤(b)中，所述关系式(四)优选按照下列表达式进行：

其中，K(k)是k时间步的卡尔曼增益，Q是系统噪声方差矩阵，q是加速度测量噪声的方差，是第k时间步的状态预测值的方差矩阵，是第k-1时间步的状态修正值的方差矩阵，H＝[1 0]表示测量矩阵，R＝r是应变测量噪声的方差。

进一步优选地，在步骤(c)中，所述关系式(五)优选按照下列表达式进行：

其中，是第k时间步状态向量的修正值。

总体而言，通过本发明所构思的以上技术方案与现有技术相比，能够取得下列有益效果：

1、本发明通过构建第k时间步悬臂梁的节点处位移x_m(k)与每个节点处的应变值与虚弯矩之间的关系式(一)，建立应变与位移之间的关系，解决了大型长细结构位移测量难度大的问题：一方面加速度和应变数据测量难度小，另一方面利用虚功原理及图乘法由应变计算得位移，不需要对大型结构进行建模获取结构振型，简单直接；

2、本发明通过构建第k时间步状态向量修正值关于状态向量预测值修正系数K(k)和节点处位移x_m(k)的关系式(五)，通过测量每个节点处的应变和加速度后即可求解获得各个节点处任意时间步的动态位移，解决了超高层结构位移测量精度低的问题，具体地，通过将加速度和应变这两种间接位移测量方法相结合，取长补短，利用加速度的高频高精度优势弥补应变的精度不足，同时利用应变修正加速度的基线偏移问题，从而提高了结构动态水平位移的精度。

3、本发明提供的方法，步骤简单，计算精度高，且抵抗测量噪声的能力强，具有较强的工程适用性。

附图说明

图1是按照本发明的优选实施例所构建的动态位移获取方法的流程图；

图2是按照本发明的优选实施例所构建的是超高层结构简化为悬臂梁模型示意图；

图3是是按照本发明的优选实施例所构建的节点角位移图；

图4是按照本发明的优选实施例所构建的虚弯矩示意图；

图5是按照本发明的优选实施例所构建的超高层结构有限元模型立面图；

图6是按照本发明的优选实施例所构建的超高层结构有限元模型截面图；

图7是按照本发明的优选实施例所构建的两种方法的位移计算结果比较图；

图8是按照本发明的优选实施例所构建的3个不同楼层两种方法的位移计算误差比较图；

图9是按照本发明的优选实施例所构建的噪声20％情况下两种方法位移计算误差比较图。

具体实施方式

为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。此外，下面所描述的本发明各个实施方式中所涉及到的技术特征只要彼此之间未构成冲突就可以相互组合。

图1是按照本发明的优选实施例所构建的动态位移获取方法的流程图，如图1所示，一种基于卡尔曼滤波数据融合的动态位移获取方法，该方法包括下列步骤：

(1)基于结构整体变形与局部应变间的关系，利用虚功原理及图乘法，由结构沿长度方向的应变计算得到垂直于长度方向的结构位移：

(1.1)：图2是按照本发明的优选实施例所构建的是超高层结构简化为悬臂梁模型示意图，如图2所示，将细长结构简化为一个竖向悬臂梁，将悬臂梁划分为n个子单元，共n+1个节点(编号0,1,…,n)；

(1.2)：在结构n+1个节点的两侧分别安装应变传感器，测量各节点的动态应变响应；

(1.3)：图3是是按照本发明的优选实施例所构建的节点角位移图，如图3所示，根据节点两侧应变差值，计算节点的角位移：

式中，s为节点处结构的宽度，Δε(k)表示节点第k时间步两侧应变的差值，dh为节点的高度，dθ(k)为节点在第k时间步的角位移；

(1.4)：图4是按照本发明的优选实施例所构建的虚弯矩示意图，如图4所示，基于虚功原理，得到悬臂梁结构各节点的虚弯矩值：

式中，是第i个节点的虚弯矩，h_j表示第j个子单元的长度；

(1.5)：根据子单元上下两端节点的角位移和虚弯矩，利用图乘法计算子单元上、下节点之间的相对位移：

式中，Δx_i为第i个子单元在第k时间步上、下节点之间的相对位移；

(1.6)：将各子单元的相对位移累加，得到结构第i个节点的水平位移：

式中，x_m(k)表示结构顶部节点第k时间步的水平位移，p是第i个节点前节点的总数量。

(2)根据位移与加速度之间的关系建立结构位移状态空间方程：

(2.1)：建立第k时间步位移与加速度之间的关系式，以及速度与加速度之间的关系式：

式中，x(k)为节点在第k时间步的位移，分别为节点在第k-1时间步的速度和加速度，Δt为两个时间步间的时间间隔，w_d为位移系统误差，w_v为速度系统误差；

(2.2)：建立真实位移与实测位移之间的关系式：

x_m(k)＝x(k)+v(k)

式中x_m(k)为位移测量值，x(k)为位移真实值，v(k)表示测量噪声；

(2.3)：将位移和速度作为状态变量，引入状态向量

式中，X(k)表示状态向量，上角标T表示向量的转置；

(2.4)：基于状态向量，建立位移的状态方程和测量方程：

x_m(k)＝HX(k)+v(k)

式中，表示状态转移矩阵，表示控制项输入矩阵，H＝[1 0]表示测量矩阵，表示系统噪声矩阵。

(3)基于状态空间方程，通过卡尔曼滤波将由应变计算所得位移与加速度进行融合，得到精度更高的结构位移：

(3.1)：设定状态向量的初值包括位移初值和速度初值；设定状态向量的预测值与修正值的方差初值

式中，和分别表示状态向量的预测值和修正值，和分别为状态预测值和修正值的方差矩阵；

(3.2)：基于位移状态方程，通过第k-1时间步的状态向量和加速度测量值预测第k时间步的状态向量：

式中，k＝2，3，…，N表示第k步迭代，为加速度测量值；

(3.3)：计算状态向量预测值的方差，用于计算卡尔曼增益，具体如下：

式中，和分别为状态预测值和修正值的方差矩阵，系统噪声方差矩阵

(3.4)：计算卡尔曼增益，作为对状态向量预测值进行修正的修正系数，具体如下：

式中，K为卡尔曼增益，测量噪声方差R＝E[vv^T]＝r；

(3.5)：基于应变计算得到的位移值对状态向量的预测值进行修正得到状态向量的修正值，具体如下：

式中，表示状态向量的修正值，x_m为由应变计算得到的位移值；

(3.6)：计算状态向量修正值的方差矩阵，具体如下：

(3.7)：步骤(3.2)至步骤(3.6)为一次卡尔曼滤波迭代过程，重复所述迭代过程，得到每一时间步的状态向量修正值，所述状态向量修正值中的位移值即为最优位移值。

图5是按照本发明的优选实施例所构建的超高层结构有限元模型立面图，图6是按照本发明的优选实施例所构建的超高层结构有限元模型截面图如图5和6所示：以超高层有限元模型为对象，来描述基于卡尔曼滤波的结构动态位移计算方法过程。超高层有限元模型为框架-核心筒结构，结构高度288m，共64层，每层层高4.5m，选择第1、7、14、24、28、35、42、49、56、64层这10个楼层作为应变测量楼层。结构外框架尺寸为50×50m，核心筒尺寸为30×30m。结构激励方式为在其顶部X方向设置水平初始位移，通过动态分析得到结构的应变、加速度和位移响应。将结构水平位移计算值作为位移基准值。

第一步：利用虚功原理及图乘法，由结构竖向应变计算得到结构水平位移。提取所述10个楼层外框架两边的Z方向应变响应，用于计算结构第64层的动态位移。

第二步：通过卡尔曼滤波将由应变所得位移与加速度进行融合。提取结构第64层的加速度，用于与应变计算所得位移进行融合。将应变所得位移、数据融合所得位移与位移基准值进行比较。

为验证本发明适用于计算超高层结构各个高度的水平位移，另选取第28层和第42层进行位移计算并将结果与基准值进行比较。图7是按照本发明的优选实施例所构建的两种方法的位移计算结果比较图，图7给出了第64层的仅用应变计算所得位移与本发明所得位移计算值与基准值之间的对比，图8是按照本发明的优选实施例所构建的3个不同楼层两种方法的位移计算误差比较图，图8给出了第64、42、28层两种位移计算值与基准值之间误差的比较结果。由图可知，本发明提出的位移计算方法误差显著小于仅使用应变数据的计算方法，这说明了本发明能有效提高位移精度。

为验证本发明所述方法的抗噪性能，对应变和加速度数据加入四种不同噪声水平的高斯白噪声模拟实际测量值，分别为5％、10％、15％、20％。拟加入的噪声数据模型为：

测量值＝awgn(准确值，SNR)

式中，awgn()为MATLAB软件中给数据加入高斯白噪声的程序，SNR为真实信号和噪声信号的比值，单位为dB。5％、10％、15％、20％高斯白噪声所对应的SNR值分别为26、20、16、14。

本发明所述位移计算误差为位移计算值与位移基准值之间的相对差异，具体如下：

式中，x_cal(k)为各时间步位移计算值，x_ref(k)为各时间步位移基准值。

表1给出了不同噪声水平下仅用应变数据的位移计算方法的误差，表2给出了不同噪声水平下本发明方法的位移计算误差。由表1可知，噪声水平从5％增大到20％过程中，仅用应变数据的位移计算方法误差从2.27％增大至8.71％。由表2可知，噪声水平从5％增大到20％过程中，本发明所述方法的误差仅从0.66％增大到1.60％。图9是按照本发明的优选实施例所构建的噪声20％情况下两种方法位移计算误差比较图，图9给出了噪声20％情况下两种方法位移计算误差比较图，图中显示卡尔曼滤波计算位移结果和基准值基本吻合，而应变方法的位移计算误差较大。由图表可知，本发明所述方法具有良好的抗噪性，在高噪声水平下本发明仍能保持高精度。

表1不同噪声水平下仅用应变数据的位移计算误差

表2不同噪声水平下本发明的位移计算误差

本领域的技术人员容易理解，以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围之内。