CN111504278A - 基于自适应频域积分的海浪检测方法 - Google Patents

Abstract

本发明的目的在于提供基于自适应频域积分的海浪检测方法，针对传统测量装置中垂直运动的被积分加速度不准确，积分位移精度低、可靠性差、易发散等问题。提出了AHRS，将载体坐标系加速度转换到地理坐标系下的天向加速度，以获得更准确的被积分加速度。AHRS还实现了对磁力计的实时在线校正，可以获得更为准确的航向角，从而提高波向的测量精度。本发明可以根据海浪运动情况自动调节滤波截止频率，提高滤波有效性及测量精度。

Description

基于自适应频域积分的海浪检测方法

技术领域

本发明涉及的是一种海洋环境监测方法，具体地说是海浪运动参数波高、周期、波向的监测方法。

背景技术

随着地球上人口的不断增长，加剧了资源的过渡开采，资源短缺已经成为了世界性难题。

而海洋中蕴含丰富的石油、天然气等自然资源，但目前人类所能使用的海洋资源十分有限。人类如果想要开发海洋资源，就需要认识海洋、对海洋进行监测，掌握海浪的相关参数。所以对海浪的检测具有重要价值和意义。

海浪检测装置是一种可以测量海浪波高、周期、波向的设备，具有灵活的工作方式，可以自动、单次、定时、定点测量海浪的运动情况。根据测量原理可以将海浪检测装置分为：声学类装置、雷达遥感类装置、加速度传感器类装置、GPS类装置，不同装置各有优缺点。综合考虑测量精度、稳定性、响应速度等性能指标，现阶段，只有加速度传感器类和GPS类装置被广泛应用于海浪的测量。虽然GPS类装置测量精度相对较高，但GPS信号易受天气影响，并且造价相对较高，不便于区域组网监测海浪情况，在测量精度要求不高的情形下，大多采用加速度传感器类的检测装置。

加速度传感器类的装置，主要通过三轴加速度传感器，测得装置三轴加速度，可以很好的反映海浪三维运动情况。将加速度信号进行二次积分转换为位移信号，从而解算出周期、波高。并通过电子罗盘计算出系统的地理航向，根据载体的水平运动情况计算出海浪的实际波向。加速度类传感器类检测装置测量的关键与难点在于加速度信号的利用与转换。传感器直接获得的是载体坐标系下的加速度，而海浪运动的周期与波高是基于地理坐标下的海浪垂直运动的位移而定义，随着海浪的运动，载体坐标系与地理坐标系往往存在偏差，传感器直接测量的加速度无法有效的表征海浪真实运动情况。加速度信号积分成位移信号的方法大致可以分为两类：时域积分、频域积分。传感器测量的加速度信号中往往掺杂着各类干扰信号，受直流信号的影响，时域积分方法在时域内的积分会造成误差的不断积累，导致积分的位移不准确，需要使用复杂的误差补偿方法进行修正。而直流等干扰在频域积分内被视为常量，不会引起误差的积累。在频域积分时受各类干扰信号影响，需要将干扰信号进行滤波消除，才可获得时域内的位移信号。但由于海浪运动周期的不确定性，固定参数的滤波器往往精度较低，甚至引起滤波异常，造成海浪的位移信号的缺失，降低了测量精度。如何设计有效的滤波器，自适应海浪周期的变化，是提高频域积分有效性的关键。并且受载体航向的测量精度影响，航向的测量易受外界干扰，导致测量精度不高。

发明内容

本发明的目的在于提供可以根据海浪运动情况自动调节滤波截止频率，提高滤波有效性及测量精度，并且通过对磁力计的在线补偿，提高载体地理航向测量精度，进而提高波向解算精度的基于自适应频域积分的海浪检测方法。

本发明的目的是这样实现的：

本发明基于自适应频域积分的海浪检测方法，其特征是：

(1)降频采集并存储MEMS输出的陀螺仪、加速度计、磁场数据；

(2)利用加速度计及磁场数据通过AHRS解算出载体姿态并进行存储；

(3)使用姿态转换矩阵将载体坐标系下加速度信号，转换到地理坐标系下，其中天向加速度即为海浪运动的垂直加速度；

(4)使用快速傅立叶变换将时域内的海浪运动的垂直加速度信号转换到频域内；

(5)将频域内的加速度信号进行二次积分成位移信号，使用自适应频域积分算法，确定上、下限截止频率，滤除相关干扰信号；

(6)解算出海浪运动的特征值：包括最大波高、最大周期、十分一波高、十分之一周期、三分之一波高、三分之一周期、平均波高、平均周期；

(7)根据地理坐标下的水平面合加速度方向及AHRS解算出的航向角，获取海浪运动的主波向和实际运动波向；

(8)将海浪的测量时间、采样频率、有效波个数、波高、周期、波向进行存储，通过通讯接口发送海浪的特征值。

本发明还可以包括：

1、将载体坐标系加速度信号转换到地理坐标系，具体为：

MEMS直接测量的数据为载体坐标下的加速度信号，并不能直接用于位移信号的直接积分，其中载体坐标系b定义为右-前-上，地理坐标系n的定义为东-北-天，

其中

为b系到n系的姿态转换矩阵，s、c代表正余弦函数，C表示转动的方向余弦矩阵，定义载体坐标系、地理坐标系下的加速度分别如下：

xyz表示载体坐标系下右前上分量，enu表示地理坐标系下的东北天分量，地理坐标系下的加速度为：

其中天向加速度

即为海浪运动的垂直加速度信号，二次积分为位移信号。

2、所述的自适应频域积分算法，根据不同的海浪运动情况，实时在线估计出上、下限截止频率：

获取加速度信号x(n)的频域变换X(k)，以频率f为横坐标，X(f)为纵坐标，得到幅值谱，对于X(kΔf)，取X(k)最大值时，对应的频率即kΔf＝f₀，即为海浪周期运动的频率f₀，频谱的横坐标取值范围为f∈[0,Fs/2]，传统的固定参数的滤波器如下：

设计的自适应滤波器如下：

其中，β为自适应因子，从而得到海浪运动的垂直位移。

本发明的优势在于：本发明可以根据海浪运动情况自动调节滤波截止频率，提高滤波有效性及测量精度。并且通过对磁力计的在线补偿，提高载体地理航向测量精度，进而提高波向解算精度。

附图说明

图1为本发明海浪检测系统的结构框图；

图2为本发明的流程图；

图3为本发明MEMS加速度输出数据频谱图；

图4为本发明水平加速度坐标系；

图5为本发明波向分布图；

图6为本发明时域积分与频域积分误差对比图；

图7为本发明姿态转换矩阵换与三轴合成加速度频域积分误差对比图；

图8为本发明频域积分自适应滤波与固定参数滤波误差对比图；

图9为本发明波高、周期测量误差图表；

图10为本发明波向测量误差图表。

具体实施方式

下面结合附图举例对本发明做更详细地描述：

结合图1-10，本发明基于自适应频域积分的海浪检测方法，其具体步骤如下：1.降频采集数据

由海浪运动测量的统计特性可知，海浪特征值的解算需要采集一段时间内的传感器信息，然而，传感器的采样频率一般都较高，如果长时间采集数据，将对处理器产生计算、处理负担，因此降频采集数据，可以在满足测量精度的同时，降低计算量，提高系统的实时性。

2.AHRS计算载体姿态

通常描述载体运动的一组欧拉角也称为姿态角包括：航向角、俯仰角、横滚角。其中航向角ψ定义北偏东为正，角度范围0～360°，俯仰角θ定义为抬头为正，角度范围-90°～90°横滚角γ定义为当运载体向右倾斜时角度为正，角度范围-180°～180°，本发明设计了AHRS系统用于实时解算载体姿态。

(1)俯仰角、横滚角解算

在AHRS系统中主要用四元数来描述刚体转动或姿态变换，四元数微分方程的矩阵表达方式如下：

其中，q₀,q₁,q₂,q₃为四元数的组成元素，在MEMS系统中

根据四元数与姿态角的关系可以得：

MEMS的陀螺漂移较大，单独工作时会逐渐发散，为了解决这个问题，引入姿态反馈作为外部量测量进行组合，以减小姿态发散程度。其中，量测姿态分别选取下文磁力计解算的航向角ψ和重力匹配解算的俯仰角θ横滚γ公式如下：

其中，

为载体坐标下的加速度三轴分量。状态量选取为：

X＝[q₀ q₁ q₂ q₃ ε_x ε_y ε_z]^T (4)

上式中ε_x,ε_y,ε_z为陀螺的常值漂移，补偿后的四元数微分方程为：

其中陀螺常值漂服从一阶马尔科夫分布，β为反相关时间常数，w(t)为激励

高斯白噪声。

因此，可以得到系统的状态方程为：

其中F(t)，G(t)为与时间相关的常数，w(t)高斯白噪声向量。F为：

系统状态方程离散化为：

X_k＝Φ_k/k-1X_k-1+w_k-1 (8)

其中Φ_k/k-1为一步状态转换矩阵，w_k-1为系统噪声。

系统的量程方程为

Z_k＝H_kX_k+v_k (9)

根据四元数与姿态角的关系，可以得到姿态角量测量对应的四元数为：

选取四元数作为量测状态量z，观测矩阵H分别为：

假设系统噪声与量测噪声满足如下关系：

因此，X的卡尔曼滤波方程如下：

(2)航向角解算

本发明采用MEMS三轴磁力计进行航向解算，对于水平静止的载体航向角计算相对比较容易公式为：

其中，B_x和B_y为校正前磁力计X轴，Y轴输出。对于运动的载体来说航向角的计算需要俯仰角和横滚角进行补偿，B_xhB_yh为补偿后磁力计输出，具体公式如下：

(3)磁力计在线校正

由于传感器自身工艺及外界磁场干扰等限制，磁力计需要经过校准才可以进行使用，其误差主要包括：硬磁/软磁偏差、非正交误差、比例误差与零偏误差、量测噪声。常用的磁力计测量模型为：

其中B_r,k为三轴磁力计校正前k时刻输出矢量，S_M＝diag(s_x,s_y,s_z)为比例误差，C_NO为载体坐标系的非正交误差阵。C′_SI＝I+C_SI，I∈i³为单位矩阵，C_SI为软磁变换系数矩阵，一般假设其为对称矩阵。

为惯性坐标系到载体坐标系的旋转矩阵，b′_HI＝C′_SIb_HI,n′_e,k＝C′_SIn_e,k，b_HI为载体自身硬磁材料产生的磁场。b_M为固定零偏矢量，n_e,k,n_s,k为k时刻外部环境噪声和传感器自身噪声。

为k时刻惯性坐标系下地磁场矢量。校正后的磁力计测量模型为：

其中，C＝S_MC_NOC′_SI，总偏移矢量b＝S_MC_NOb′_HI+b_M，测量误差n_k＝S_MC_NOn′_e,k+n_s,k。S_M,C_NO,C_SI,b_HI,b_M等参数不需要，也不能被显示求解，因此，求得T和b

根据测量点的地磁场模值B_m为常量的原则，可以推导出校正参数状态空间模型：

令T^TT＝G代入上式可得，

则对{T,b}的估计转换为对{G,b}的估计。其中G为对称矩阵，假设其参数为

设b＝[b_x b_y b_z]^T，状态变量x＝[A B C D E F b_x b_y b_z]^T可得状态方程为

x_k＝x_k-1 (19)

根据式(18)构造状态变量x的隐式量测方程：

其中，

h₁(x_k)＝L_kx_k (21)

为校正前磁力计三轴输出量。v_k为噪声项，n_k～N(0,N_k)，

为测量噪声协方差。由于式(20)为非线性方程，因此需要使用非线性滤波器对对{G,b}进行估计。本发明采用一种实时性较好的非线性滤波器EKF。

在

和n_k＝0对量测方程进行一阶泰勒展开得

可以得到线性量测方程：

分别为h对x_k和n_k在

n_k＝0处的雅可比(Jacobian)矩阵，并且有

M_kn_k～N(0,r_k) (29)

其中

因此x_k的EKF滤波公式如下：

随着磁场量测值

进行迭代更新，当时间平均残差

收敛时，迭代更新结束即：

其中ε为小正数，

当残差收敛时，得到的估计值可以用于下一步的校正计算。由此可以得到{G,b}对矩阵G进行特征值分解可得

G＝QΣQ^-1 (32)

其中，Σ＝diag(λ₁,λ₂,λ₃)为对角阵，λ₁λ₂λ₃为G矩阵对应的特征值，Q＝[α₁,α₂,α₃]，又因为G＝T^TT由此可以解得

所以可以获得校正参数{T,b}利用式(5)即可获得校正后的磁场

再根据式(14)与式(13)即可解算出载体的航向角ψ。

3.将载体坐标系加速度转换到地理坐标系

MEMS直接测量的数据为载体坐标下的加速度信号，并不能直接用于位移信号的直接积分，其中载体坐标系b定义为右-前-上，地理坐标系n的定义为东-北-天是常用的导航坐标系。

其中

为b系到n系的姿态转换矩阵，s,c代表正余弦函数，C表示转动的方向余弦矩阵。定义载体坐标系，地理坐标系下的加速度分别如下：

xyz表示载体坐标系下右前上分量，enu表示地理坐标系下的东北天分量。因此，地理坐标系下的加速度为：

其中天向加速度

即为海浪运动的垂直加速度信号，二次积分为位移信号。

4.频域积分

频域积分和时域积分的主要区别在于频域上正余弦信号的关系，频域内积分本质上就是正弦余弦信号的相互转换从而简化积分计算。因此，可以有效地减少干扰信号在时域积分过程中带来的误差积累影响。具体积分步骤如下：

(1)使用快速傅立叶变换(FFT)将时域信号转换到频域，保证精度的同时可

以极大提高计算效率，FFT的公式如下：

假设时域内振动信号的离散数据{x(n)}(n＝0,1,2,……,N)令

x(n)＝x₁(n)+x₂(n) (38)

x₁(n)和x₂(n)的长度都是N/2，x₁(n)是偶数列，x₂(n)是奇数列。

X₁和X₂分别为x₁(n)和x₂(n)对应的傅立叶变换(DFT)。

(2)频域内一次积分计算公式为

(3)频域内二次积分计算公式为

其中

式中X(k)为x(n)的FFT，Δf为频率分辨率，H(k)为带通滤波器滤除相关干扰信号，f_min和f_max为下限和上限截止频率，y(n)为时域内的位移信号，加速度信号在频域内二次积分后过带通滤波器后，根据傅立叶逆变换公式将频域内的位移信号转换为时域内的位移信号y(n)。5.频域积分自适应滤波算法

采用加速度传感器类测量海浪特征值(波高、周期、波向)时，由于加速度计的噪声特性及实际测量环境中存在的噪声等误差，会使海浪特征值测量产生一定的偏差，所以需要数字滤波器进行滤波处理。时域积分时，传统的数字滤波器大多都采用固定的上下限截止频率，常用的带通滤波器如式(42)所示：

这些截止频率往往与实际存在一定偏差，尤其在频域积分内，影响尤为显著。频域积分时由于海浪运动周期范围为2-30s，不同周期的海浪运动对应的上下限截止频率大不相同，如果采用固定的上下限截止频率，滤波效果将大大降低，甚至滤波发散。

(1)加速度信号频谱图

根据频域积分的特点，可以获得加速度信号x(n)的频域变换X(k)。如图3所示，为某实测加速度信号频谱图，以频率f为横坐标，X(f)为纵坐标，即可得到幅值谱。由采样定理可知，可以正确采集信号的频率不超过Fs/2，Fs为采样频率。FFT中频率超过Fs/2的部分是不准确的，因此频谱的横坐标取值范围为f∈[0,Fs/2]。

(2)估计海浪运动的频率

对于X(kΔf)，取X(k)取最大值时，即海浪运动的主要频率范围，对应的频率即kΔf＝f₀，即为海浪周期运动的频率f₀，其最大误差为Δf/2。即使加速度信号中混有噪声，噪声也会分布在频谱的不同频率上，并且幅值低于海浪主要运动的峰值，因此上述方法仍可估计海浪运动的真实频率。

(3)估计滤波器上下限截止频率

如图3频谱图所示，加速度计的采样频率为4Hz，所以频谱图的横坐标范围为0～2Hz。由海浪的相关理论可知，海浪的周期一般为10s，因此，海浪的频率范围在0.1Hz左右。根据局部放大图可知，当加速度计频率小于0.1Hz的信号为加速度计的零偏误差，海浪实际测量过程中的误差主要集中在低频段，所以加速度计频率大于0.1～2Hz的信号为海浪的有用信号，因此滤波器下限截止频率f_min取值尤为重要。但实际海浪测量时，一段时间内海浪的周期并不是恒定不变的，因此其频谱图是在幅值最大值对应的频率f₀附近波动的，如果f_min取值偏大，会造成积分出的位移信号缺失，f_min取值偏小，位移信号会包含加速度计零偏等低频误差，考虑以上因素，设计了如下自适应滤波器：

其中，β为自适应因子，由实际工程经验获得。该自适应滤波器具有很好的实时性及有效性，可以有效的估计出f_min与f_max，提高滤波精度与可靠性，从而得到海浪运动相对准确的垂直位移。

6.波高、周期解算

《海滨观测规范》可得，海浪的统计方法以水平位移的平均值为零线，由零点定理可知，存在海浪运动的两个相邻的位移点y₁(n₁)和y₂(n₂)使得y₁y₂<0，则任取n₁、n₂为零点对应的时刻，其误差小于等于Fs/2。采用上跨零点法，定义海浪上升与零点交点作为起始点，当海浪运动到零点以下，再次上升与与零点的交点作为终点。起始点与终点之间的时间差为周期，起点与终点之间的波峰与波谷的垂直高度差为波高。记录一个观测周期内n个有效波的波高和周期，并对波高进行从大到小排序：H₁,H₂,H₃......H_n，及其对应周期为T₁,T₂,T₃......T_n。

海浪运动的波高和周期特征值包括：最大波高H_max，及其对应的周期T_max、十分之一波高H_1/10及其对应的周期T_1/10、有效波高H_1/3及其对应的周期T_1/3、平均波高H_mean，及其对应的周期T_mean。其计算公式如下：

7.波向解算

将一个海浪周期内跨越零点时地理坐标系下的水平合成加速度作为海浪的运动方向α。以加速度载体坐标系Ay轴作为正北方向，以Ax轴作为正东方向建立如图4所示坐标系，其中水平合加速度与y轴的夹角为α即与北向夹角，由前文可知，AHRS可以获得载体地理坐标系下实际地理航向角ψ，根据地理坐标系下水平合加速度与地理航向关系即可得到海浪运动的地理波向σ₀，即水平合运动方向与北向的夹角。具体计算公式如下：

其中

和

为零点处地理坐标系下的加速度分量，其中σ₀与α取值范围均为：0～360°，取北偏东为正。因此，可以得到一个周期内海浪运动的波向。由于海浪波向是一个长期统计参数，所以以22.5°为间隔将图5所示坐标系分割成16个扇区，海浪波向的出现频率最高的扇区为主波向m₁(0-15)，再将每个扇区分成22个小扇区，出现频率最高的扇区为m₂(0-21)，从而得到观测时间内的波向σ，分辨率为1°。具体公式如下：

σ＝m₁×22.5°+m₂ (46)

8.存储解算数据

为了便于对海浪运动的分析，将解算后的海浪特征值、采样时间等关键信息通过FAT系统存入到SD卡中，并通过通讯接口发送测量信息给相应的设备。

一、对本发明进行实测实验验证：

为了验证本发明海浪特征值解算精度，测试实验采用海浪测试平台，其中波高测量误差：±185um，周期测量误差：±0.05s。可以分别对1m，3m，6m三个波高点七个周期进行测量。具体测量结果如图9所示，由于工程实践中，侧重关注有效波高及周期，因此表中皆为有效波高，周期测量结果及误差。因为波浪周期波高的统计特性，在采样周期内，海浪的周期越长，速度越慢，采样的有效波个数就越少，则由于个别测量周期内误差过大，将使整体的误差变大，因此，波周期越长，测量精度相对越低，周期越短，测量精度越高。测量经多次测量，本发明波高测量精度为：0.1±5％×实际波高(m)，周期测量精度为：±0.25s。波向测量结果如表二所示，波向测量精度为：±5°，分辨率为1°。

图6为时域积分与频域积分误差对比图，相对时域积分，频域积分测量结果更为准确。频域内的误差量为常值，并不会随着积分次数的增加而逐渐积累，而时域积分，逐渐呈现发散趋势，充分说明了频域积分方法的有效性。

图7为加速度通过坐标转换与三轴合成加速度的误差对比图，由图可知，坐标转换后的加速度积分的波高误差相对更小，说明其更真实的反映了垂直方向的运动情况。

图8为自适应滤波与传统频域积分滤波误差对比图，可以看出频域积分自适应滤波算法，具有更高的滤波精度、收敛速度快、稳定性好，充分说明其估计的滤波器上下限截止频率，更加符合实际运动情况，在滤除相关干扰误差的前提下，更好的反映运动状态。通过实测数据仿真及本装置实时输出的串口信息可知，本发明的方法测量误差在允许范围内，可以实现对海浪特征值的准确解算。

Claims (3)

1.基于自适应频域积分的海浪检测方法，其特征是：

(1)降频采集并存储MEMS输出的陀螺仪、加速度计、磁场数据；

(2)利用加速度计及磁场数据通过AHRS解算出载体姿态并进行存储；

(3)使用姿态转换矩阵将载体坐标系下加速度信号，转换到地理坐标系下，其中天向加速度即为海浪运动的垂直加速度；

(4)使用快速傅立叶变换将时域内的海浪运动的垂直加速度信号转换到频域内；

(5)将频域内的加速度信号进行二次积分成位移信号，使用自适应频域积分算法，确定上、下限截止频率，滤除相关干扰信号；

(6)解算出海浪运动的特征值：包括最大波高、最大周期、十分一波高、十分之一周期、三分之一波高、三分之一周期、平均波高、平均周期；

(7)根据地理坐标下的水平面合加速度方向及AHRS解算出的航向角，获取海浪运动的主波向和实际运动波向；

(8)将海浪的测量时间、采样频率、有效波个数、波高、周期、波向进行存储，通过通讯接口发送海浪的特征值。

2.根据权利要求1所述的基于自适应频域积分的海浪检测方法，其特征是：将载体坐标系加速度信号转换到地理坐标系，具体为：

MEMS直接测量的数据为载体坐标下的加速度信号，并不能直接用于位移信号的直接积分，其中载体坐标系b定义为右-前-上，地理坐标系n的定义为东-北-天，

其中

为b系到n系的姿态转换矩阵，s、c代表正余弦函数，C表示转动的方向余弦矩阵，定义载体坐标系、地理坐标系下的加速度分别如下：

xyz表示载体坐标系下右前上分量，enu表示地理坐标系下的东北天分量，地理坐标系下的加速度为：

其中天向加速度

即为海浪运动的垂直加速度信号，二次积分为位移信号。

3.根据权利要求1或2所述的基于自适应频域积分的海浪检测方法，其特征是：所述的自适应频域积分算法，根据不同的海浪运动情况，实时在线估计出上、下限截止频率：

获取加速度信号x(n)的频域变换X(k)，以频率f为横坐标，|X(f)|为纵坐标，得到幅值谱，对于X(kΔf)，取X(k)最大值时，对应的频率即kΔf＝f₀，即为海浪周期运动的频率f₀，频谱的横坐标取值范围为f∈[0,Fs/2]，传统的固定参数的滤波器如下：

设计的自适应滤波器如下：

其中，β为自适应因子，从而得到海浪运动的垂直位移。

Images