CN110533005B - 一种复杂海况下船体形变测量方法 - Google Patents

Abstract

本发明公开了一种复杂海况下船体形变测量方法，包括步骤：数据采集；基于汉宁窗的FIR低通滤波器设计；状态向量降维处理；状态方程和量测方程建立；卡尔曼滤波估计。本发明可以解决复杂海况下的船体动态形变和静态形变估计问题。相比于目前普遍采用的基于角速度匹配的方法，解决了计算量大，估计精度低的问题，精度和收敛性均得到进一步提高。且在复杂海况下能够降低外界噪声干扰，提高形变估计的稳定性。

Description

一种复杂海况下船体形变测量方法

技术领域

本发明属于舰船检测技术，具体涉及一种复杂环境下船体形变测量方法。

背景技术

随着海军装备的不断发展和作战需求的提升，各种舰载武器也朝着智能化、高精度、高水平方向发展。目前舰船一般装有各种攻击和防御设备以及导航***等，比如跟踪***、导弹发射***、雷达***、瞄准***等。不同的舰载设备又根据作战任务不同分布于舰船不同战位点，为了保证这些设备的统一协调运行，舰船***需向这些设备提供高精度的导航信息如：位置、速度和姿态信息，以保证建立起全舰统一的空间基准。然而，船体本身并非一个绝对刚体，由于受到海浪等冲击导致这种非刚体特性使舰船在不同海况下产生三轴方向的形变角。形变角导致舰船设备相对于舰船摇摆中心产生基准误差。因此，舰船上的各个设备***也就不能得到精确的导航参数。这就降低了武器的命中率。所以为保证各个子***间协调工作时的精度，精确测量出船体变形是非常有必要的。

目前船体形变测量最常用的方法是角速度匹配测量法，其工作原理为：在舰船上两个不同的部位安装两个惯性测量***(IMU)，当这两个部位之间没有变形角时，坐标系是吻合的，两个IMU的输出是一致的。但是若两个部分之间存在变形角时，两个坐标系将不再吻合，此时两个IMU所敏感的载体运动如纵摇、横摇和航向也不再一样。利用两个IMU输出的陀螺角速度之差当做观测量，就可以采用Kalman滤波等算法来估计船体变形角。

针对光纤陀螺船体形变测量***在复杂海况下测量精度低的问题，提出了一种改进的角速度匹配船体形变测量方法。该方法在传统的基于角速度差作为观测量的方法基础上，通过优化***模型，降低维数，提出了一种简化模型，并建立线性量测方程。在此基础上，考虑在复杂海况下基于角速度的惯性匹配船体形变测量方法受到环境干扰较强，此时测量精度将降低，主要的原因如舰船受海浪引起的高频振颤等致使形变测量方法滤波过程中受到有色噪声干扰，加之存在***建模误差，从而导致滤波器性能降低，严重时可能导致滤波器发散。引入FIR低通滤波器，保证复杂海况下船体形变估计精度。

发明内容

发明目的：本发明的目的在于解决复杂海况下船体形变测量易受噪声干扰，精度低的缺点，提供一种复杂海况下船体形变测量方法。

上述的目的通过以下技术方案实现：

一种复杂海况下船体形变测量方法，该方法包括如下步骤：、

步骤1：在船体的中心位置安装主IMU，在船体的边缘位置安装子IMU，主IMU用于采集光纤陀螺信号

包括

其中，

表示主IMU测量的x轴角速度值；

表示主IMU测量的y轴角速度值；

表示主IMU测量的z轴角速度值；子IMU用于采集的光纤陀螺信号

包括

其中，

表示子IMU测量的x轴角速度值；

表示子IMU测量的y轴角速度值；

表示子IMU测量的z轴角速度值；***根据主IMU和子IMU采集到的信号，输出船体动态形变角(λ_x,λ_y,λ_z)；船体静态形变角

其中，λ_x表示船体在x轴方向的动态形变角；λ_y表示表示船体在y轴方向的动态形变角；λ_z表示船体在z轴方向的动态形变角。

表示船体在x轴方向的静态形变角；

表示表示船体在y轴方向的静态形变角；

表示船体在z轴方向的静态形变角；

步骤2：基于汉宁窗的FIR低通滤波器设计：

汉宁窗函数数学表达式为：

式中：W(n)表示输出序列；n表示序列元素个数；N表示窗口值。

步骤3：子IMU光纤陀螺数据低通滤波处理：对子IMU的三个轴的陀螺输出同时进行基于汉宁窗的FIR低通滤波器处理；

步骤4：状态向量降维处理及建立状态方程；

步骤5：建立量测方程

步骤6：卡尔曼滤波器设计。

所述的复杂海况下船体形变测量方法，步骤4中所述的状态向量降维处理及建立状态方程的具体步骤是：

将主IMU测量数据的坐标系定义为m坐标系；将子IMU测量数据的坐标系定义为s坐标系，在静态情况下，主IMU和子IMU之间存在失准角为

该失准角为静态形变角；当船受到影响时，主IMU和子IMU之间的形变角度发生改变为动态形变角，记为λ，且动态形变角速率为μ，则有

则两处测量点的形变角为：

根据船体形变角产生的原因，将船体形变角分为静态形变角和动态形变角两部分；其中，静态形变角，表示为：

其中，x表示三维空间中x轴方向；y表示三维空间中y轴方向。由于静态形变角在一定时间内可以认为是一个固定不变的常量，因此进一步静态形变角速度表示为：

动态形变用二阶马尔卡夫过程作为船体变形的数学模型，动态形变角自相关函数表示为：

其中i＝x,y,z，μ_i为不规则系数表征了动态形变角的复杂程度，χ_i是动态变形的主频率，D_i为动态形变角的方差；上述相关函数对应的变形滤波器方程表示为：

其中，w(t)为白噪声，其均值为0，方差为1；

根据统计规律，陀螺漂移由常值漂移和随机漂移两部分构成，其中，常值漂移被认为是为一个常数，即：

ε₀＝0

随机漂移ε_r被认为是受白噪声w(t)驱动的一阶马尔科夫过程，表示为：

其中，τ_r是马尔可夫过程的相关时间，w_ε是马尔可夫过程白噪声；

建立状态方程为：

其中，X表示状态向量；F表示***状态转移矩阵；G表示***噪声矩阵；W表示***噪声向量。传统的基于角速度匹配的船体形变测量方法中状态向量X为：

其中，ε_1c＝[ε_1cx ε_1cy ε_1cz]^T为主IMU陀螺常值漂移；ε_1r＝[ε_1rx ε_1ry ε_1rz]^T为主IMU陀螺随机漂移；ε_2c＝[ε_2cx ε_2cy ε_2cz]^T为子IMU陀螺常值漂移；ε_2r＝[ε_2rx ε_2ry ε_2rz]^T为子IMU陀螺随机漂移。

考虑到主IMU陀螺误差参数和子IMU陀螺误差参数不可观的因素，以及进一步降低该方法的计算量，将状态向量X降低到9维：

其中，

表示船体静态静态形变角；λ表示船体动态形变角；μ表示船体动态形变角速度；则对应的F***状态转移矩阵；G***噪声矩阵表示如下：

所述的复杂海况下船体形变测量方法，步骤5中所述的建立量测方程的具体步骤如下：

量测方程考察主IMU和子IMU的陀螺输出角速度之差作为观测量，设主IMU和子IMU的输出分别为

和

由于主IMU和子IMU安装位置存在形变角，因此用姿态转移矩阵

来描述两个位置的坐标系变换关系，由于形变角为小角，近似有：

其中，

和λ×分别为

和λ的反对称阵。取

为

的反对称阵，那么有：

记观测量为：

***观测方程为：

Z＝HX+V

其中，H为量测方程转移矩阵；V为量测白噪声；

降维改进后的H阵为：

所述的复杂海况下船体形变测量方法，步骤6中所述的卡尔曼滤波方程如下：

X_k,k-1＝FX_k-1

K_k＝P_k,k-1H^T(HP_k,k-1H+R)^-1

X_k＝X_k,k-1+K_k(Z_k-HX_k,k-1)

P_k＝(I-K_k)P_k,k-1。

有益效果：

与现有技术相比，本发明具有如下优点：

(1)本发明在传统的基于角速度匹配的船体形变测量方法的基础上将状态向量维数由原来的21维降低到9维，在保证精度的情况下减少了运算量。

(2)本发明将FIR低通滤波器引入到船体形变测量中，设计了适合船体形变测量的滤波器参数，基于汉宁窗的FIR低通滤波器。

(3)与传统的基于角速度匹配的船体形变测量方法相比，本发明在复杂海况下能够进一步提高船体形变测量的精度包括动态形变角和静态形变角。

附图说明

图1为本发明所描述的复杂海况下船体形变测量方法原理图；

图2为船体静态形变示意图；

图3为船体动态形变示意图；

图4为两种算法下船体静态形变角仿真结果；

图5为两种算法下船体静态形变角误差仿真结果；

图6为两种算法下船体动态形变角仿真结果；

图7为两种算法下船体动态形变角误差仿真结果。

具体实施方式

下面对本发明技术方案进行详细说明，但是本发明的保护范围不局限于所述实施例。

1、一种复杂海况下船体形变测量方法，其特征为：

已知量：主IMU采集的光纤陀螺信号

子IMU采集的光纤陀螺信号

输出量：船体动态形变角(λ_x,λ_y,λ_z)；

船体静态形变角

步骤1：采集主IMU和子IMU输出的光纤陀螺信号。

主IMU安装在船体的中心位置，这里采用高精度的光纤陀螺，可采集到三轴的光纤陀螺输出信息：

同理子IMU采集到的三轴光纤陀螺输出信息为：

步骤2：基于汉宁窗的FIR低通滤波器设计。

汉宁窗函数数学表达式为：

截至频率：船在复杂环境下受到不同风浪冲击，形变大小会出现变化，同时陀螺传感器受海浪引起的高频振颤等致使测量数据受到有色噪声干扰，从而影响船体形变的估计。而真实的陀螺数据是反映船体角运动，如船体运动和船体形变导致的陀螺数据变化，而这些都属于缓慢变化，因此陀螺的有效数据主要表现为低频特性。因此，通过以上分析，可以采用低通滤波器实现复杂环境下陀螺噪声的滤波处理。低通滤波器设计中截至频率的含义为，信号在低于该频率的频带内期望具有尽可能好的低通性能，而在高于该频率的频带内锐减。将截至频率设置为0.2Hz。

滤波器阶数：对于模型阶数的确定，根据其阶数越大精度越高，但相应计算量和时间延迟也增大的特点，结合实际滤波效果进行确定。选择100阶作为FIR低通滤波的阶数。

步骤3：子IMU光纤陀螺数据低通滤波处理。

对子IMU的三个轴的陀螺输出同时进行基于汉宁窗的FIR低通滤波器处理。

步骤4：状态向量降维处理及状态方程建立

为测量船体形变，一般在船主惯导***(MINS)(一般是主IMU组成)的安装位置为参考基准点，其测量数据的坐标系定义为m坐标系。在MINS(子IMU组成)位置可以获得该点的角速度、速度、姿态等信息。在船某形变待测处安装另一SINS，其测量数据的坐标系定义为s坐标系。在静态情况下，安装后的两处INS位置已经存在失准角为

该失准角为静态形变角。如图2所示。当船受到海浪、温度、负载等影响时，测量点的形变角度发生改变为λ，且动态形变角速率为μ，则有

如图3所示。则两处测量点的形变角为：

根据船体形变角产生的原因，可以将船体形变角分为静态形变角和动态形变角两部分。其中，静态形变角主要由于船体受热不均和长时间的结构老化所造成，其特点是周期长、变化慢，在一定时间内可以认为是一个固定不变的常量，因此，静态形变角可表示为

进一步，静态形变角速度可以表示为：

舰船在海上受到不同海况的影响，其动态形变和随机噪声驱动的随机过程相似，该过程可以用一个白噪声驱动的二阶马尔科夫过程来近似描述。所以用二阶马尔卡夫过程作为船体变形的数学模型。动态形变角自相关函数可表示为：

其中i＝x,y,z，μ_i为不规则系数表征了动态形变角的复杂程度。χ_i认为是动态变形的主要频率。D_i为动态形变角的方差。上述相关函数对应的变形滤波器方程可以表示为：

其中，w(t)为白噪声，其均值为0，方差为1。

大量的统计规律表明，陀螺漂移主要由常值漂移和随机漂移两部分构成，其中，常值漂移可认为在一段时间内不随时间发生变化，近似理解为一个常数，即：

ε₀＝0

随机漂移ε_r可以认为是受白噪声w(t)驱动的一阶马尔科夫过程，可以表示为：

其中，τ_r是马尔可夫过程的相关时间，w_ε是马尔可夫过程白噪声。

建立状态方程为：

其中，X表示状态向量；F表示***状态转移矩阵；G表示***噪声矩阵；W表示***噪声向量。传统的基于角速度匹配的船体形变测量方法中状态向量X为：

其中，ε_1c＝[ε_1cx ε_1cy ε_1cz]^T为主IMU陀螺常值漂移；ε_1r＝[ε_1rx ε_1ry ε_1rz]^T为主IMU陀螺随机漂移；ε_2c＝[ε_2cx ε_2cy ε_2cz]^T为子IMU陀螺常值漂移；ε_2r＝[ε_2rx ε_2ry ε_2rz]^T为子IMU陀螺随机漂移。

考虑到主IMU陀螺误差参数和子IMU陀螺误差参数不可观等因素，以及进一步降低该方法的计算量，将状态向量X降低到9维：

其中，

表示船体静态静态形变角；λ表示船体动态形变角；μ表示船体动态形变角速度。则对应的F***状态转移矩阵；G***噪声矩阵可表示如下：

步骤5：量测方程建立

量测方程考察两个位置IMU的陀螺输出角速度之差作为观测量。设两个IMU的输出分别为

和

由于两个IMU安装位置存在形变角，因此可以用姿态转移矩阵

来描述两个位置的坐标系变换关系。由于形变角为小角，近似有：

其中，

和λ×分别为

和λ的反对称阵。取

为

的反对称阵，那么有：

记观测量为：

***观测方程为：

Z＝HX+V

其中，H为量测方程转移矩阵；V为量测白噪声。

降维改进后的H阵为：

步骤5：卡尔曼滤波器设计。

卡尔曼滤波方程如下：

X_k,k-1＝FX_k-1

K_k＝P_k,k-1H^T(HP_k,k-1H+R)^-1

X_k＝X_k,k-1+K_k(Z_k-HX_k,k-1)

P_k＝(I-K_k)P_k,k-1

下面结合具体实施例，进一步阐明本发明，应理解这些实施例仅用于说明本发明而不用于限制本发明的范围，在阅读了本发明之后，本领域普通技术人员对本发明的各种等价形式的修改均落于本申请所附权利要求所限定的范围。

仿真参数如下：

IMU采样频率：200Hz；仿真时间：500s；

静态形变角：俯仰角0.3°，横滚角0.4°，航向角-0.2°；

动态形变角：二阶马尔科夫过程；

陀螺常值漂移：主IMU 0.005°/h，子IMU 0.01°/h；

载体的运动模型为：设载体运动的俯仰角为θ，横滚角为γ，航向角为ψ

其中，摇摆幅度依次是：θ_m＝6°，γ_m＝8°，ψ_m＝3°；初始相位依次是：θ₀＝0°，γ₀＝0°，ψ₀＝0°。

从图4到图7可以看出，本方法在复杂海况下的静态角估计精度和动态角估计精度都要优于传统的方法。

Claims (2)

1.一种复杂海况下船体形变测量方法，其特征在于，该方法包括以下步骤：

步骤1：在船体的中心位置安装主IMU，在船体的边缘位置安装子IMU，主IMU用于采集光纤陀螺信号

包括

其中，

表示主IMU测量的x轴角速度值；

表示主IMU测量的y轴角速度值；

表示主IMU测量的z轴角速度值；子IMU用于采集的光纤陀螺信号

包括

其中，

表示子IMU测量的x轴角速度值；

表示子IMU测量的y轴角速度值；

表示子IMU测量的z轴角速度值；***根据主IMU和子IMU采集到的信号，输出船体动态形变角(λ_x,λ_y,λ_z)；船体静态形变角

其中，λ_x表示船体在x轴方向的动态形变角；λ_y表示船体在y轴方向的动态形变角；λ_z表示船体在z轴方向的动态形变角；

表示船体在x轴方向的静态形变角；

表示船体在y轴方向的静态形变角；

表示船体在z轴方向的静态形变角；

步骤2：基于汉宁窗的FIR低通滤波器设计：

汉宁窗函数数学表达式为：

式中：W(n)表示输出序列；n表示序列元素个数；N表示窗口值；

步骤3：子IMU光纤陀螺数据低通滤波处理：对子IMU的三个轴的陀螺输出同时进行基于汉宁窗的FIR低通滤波器处理；

步骤4：状态向量降维处理及建立状态方程；

步骤5：建立量测方程

步骤6：卡尔曼滤波器设计；

步骤4中所述的状态向量降维处理及建立状态方程的具体步骤是：

将主IMU测量数据的坐标系定义为m坐标系；将子IMU测量数据的坐标系定义为s坐标系，在静态情况下，主IMU和子IMU之间存在失准角为

该失准角为静态形变角；当船受到影响时，主IMU和子IMU之间的形变角度发生改变为动态形变角，记为λ，且动态形变角速率为μ，则有

则两处测量点的形变角为：

根据船体形变角产生的原因，将船体形变角分为静态形变角和动态形变角两部分；其中，静态形变角表示为：

其中，x表示三维空间中x轴方向；y表示三维空间中y轴方向；z表示三维空间中z轴方向；由于静态形变角在一定时间内可以认为是一个固定不变的常量，因此进一步静态形变角速度表示为：

动态形变用二阶马尔可夫过程作为船体变形的数学模型，动态形变角自相关函数表示为：

其中i＝x,y,z，μ_i为不规则系数表征了动态形变角的复杂程度，χ_i是动态变形的主频率，D_i为动态形变角的方差；上述相关函数对应的变形滤波器方程表示为：

其中，w(t)为白噪声，其均值为0，方差为1；

根据统计规律，陀螺漂移由常值漂移和随机漂移两部分构成，其中，常值漂移被认为是为一个常数，即：

ε₀＝0

随机漂移ε_r被认为是受白噪声w(t)驱动的一阶马尔可夫过程，表示为：

其中，τ_r是马尔可夫过程的相关时间，w_ε是马尔可夫过程白噪声；

建立状态方程为：

其中，X表示状态向量；F表示***状态转移矩阵；G表示***噪声矩阵；W表示***噪声向量；传统的基于角速度匹配的船体形变测量方法中状态向量X为：

其中，ε_1c＝[ε_1cx ε_1cy ε_1cz]^T为主IMU陀螺常值漂移；ε_1r＝[ε_1rx ε_1ry ε_1rz]^T为主IMU陀螺随机漂移；ε_2c＝[ε_2cx ε_2cy ε_2cz]^T为子IMU陀螺常值漂移；ε_2r＝[ε_2rx ε_2ry ε_2rz]^T为子IMU陀螺随机漂移；

考虑到主IMU陀螺误差参数和子IMU陀螺误差参数不可观的因素，以及进一步降低该方法的计算量，将状态向量X降低到9维：

其中，

表示船体静态形变角；λ表示船体动态形变角；μ表示船体动态形变角速度；则对应的F***状态转移矩阵；G***噪声矩阵表示如下：

步骤5中所述的建立量测方程的具体步骤如下：

量测方程考察主IMU和子IMU的陀螺输出角速度之差作为观测量，设主IMU和子IMU的输出分别为

和

由于主IMU和子IMU安装位置存在形变角，因此用姿态转移矩阵

来描述两个位置的坐标系变换关系，由于形变角为小角，近似有：

其中，

和λ×分别为

和λ的反对称阵，取

为

的反对称阵，那么有：

记观测量为：

***观测方程为：

Z＝HX+V

其中，H为量测方程转移矩阵；V为量测白噪声；

降维改进后的H阵为：

2.根据权利要求1所述的复杂海况下船体形变测量方法，其特征在于，步骤6中所述的卡尔曼滤波方程如下：

X_k,k-1＝FX_k-1

K_k＝P_k,k-1H^T(HP_k,k-1H+R)^-1

X_k＝X_k,k-1+K_k(Z_k-HX_k,k-1)

P_k＝(I-K_k)P_k,k-1。

Images