陈志刚 韩金良

（1.海装沈阳局驻大连地区第三代表室 大连 116021）（2.天津航海仪器研究所 天津 300131）

1 引言

随着对导航精度和自主性要求的不断增加，地磁导航成为水下航行器自主导航的重要研究方向［1］。然而，由于地磁随位置变化幅度较小，加上人文以及工业磁干扰，导致地磁测量数据存在偏差［2～3］。因此，对地磁测量信号中的大幅值高频噪声进行有效滤波，从而进一步提高地磁匹配率和导航系统稳定性一直是地磁导航领域的研究热点之一。目前，对地磁测量信号高频降噪的方法主要有基于快速傅里叶变换的频域滤波、小波降噪和模态分解等方法［5］。小波降噪需要预先得到干扰磁场特征来确定小波基进行干扰抑制［6］。快速傅里叶变换滤波需要预先对噪声频域特征进行判断。经验模态分解具有自适应性、完备性和近似正交性的特点，通过滤除异常信号，能够有效提高导航精度［7～8］。

基于CEEMD的小波阈值去噪方法针对CEMMD分解可能造成有效信息损失以及小波阈值去噪方法存在的不足。基于CEEMD和PCA的组合去噪方法［9］，通过模拟实验取得了很好的去噪效果，并将该方法应用于探测雷达数据中，验证了该方法的可靠性。基于压缩感知理论小波阈值去噪方法，可以在保证高信噪比的基础上更为有效地保留弱有效信号［10］。多尺度多方向主成分分析去噪的CEEMD方法对高频GNSS同震位移进行去噪，实验结果表明该方法可以有效地削弱低频系统误差和高频白噪声，提高GNSS定位的精度［11］。由于地磁测量信号多存在高幅值的高频干扰，使用CEEMD分解会产生模态混叠现象［12］。形态滤波方法是基于积分几何、随机集合论等数学理论建立起来的一种非线性信号处理方法，该方法能有效抑制大地电磁信号中的强噪声干扰［13～15］。本文结合形态滤波与CEEMD，使用CEEMD对预处理信号进行模态分解，针对IMF的叠加选择问题，提出一种基于归一化自相关函数的模态选择方法。

2 基本原理

形态滤波方法的主要思想包括形态变换运算方式和结构元素选取，以一维离散信号为例，设f(n)，n=1，2，...，N为长度为N的原始测量信号，g(m)，m=1，2，...，M为长度为M的结构元素，且N＞＞M。则关于f(n)和g(m)的腐蚀变换为

从式（1）和式（2）中可以看出，腐蚀相当于求离散函数在滑动滤波窗即结构元素内的极小值，用来剔除边界不平滑的凸起部分，使目标收缩，孔洞扩张；膨胀运算则与腐蚀运算相反，相当于求f(n)在结构元素内的极大值，加大了信号的估值，扩展了峰顶。开、闭运算是最基本的形态滤波器，是在腐蚀和膨胀的基础上衍生而来的复合运算，其中，开运算(°)为

开运算对原信号腐蚀后膨胀，以消除原信号中的尖峰细节，使信号轮廓光滑，抑制正脉冲噪声。闭运算对原信号先膨胀后腐蚀，以弥补信号谷底，滤除低谷噪声。在实际应用中，往往结合二者共同用来抑制噪声。

采用同一结构元素，通过级联开、闭运算，定义了形态开-闭（OC）和闭-开（CO）滤波器：

为减少输出上的畸变，常使用OC和CO组合对信号进行处理：

式中，y(n)为组合形态滤波器的输出结果。

广义形态OC和CO滤波器的数学描述分别为：

则广义形态滤波器输出为

通过广义形态滤波器滤除地磁测量噪声的强干扰脉冲，得到预处理信号，使用CEEMD分解进一步将噪声分量与地磁信号分量分离，得到稳定的地磁场信息。

CEEMD是在EMD分解的基础上，通过对原信号多次叠加正负对应的白噪声，补充原始信号的频率尺度，减少包络线的拟合误差。具体步骤如下：

1）对原信号第i次叠加正负高斯白噪声，得到两组模态分量：

其中，X为原始信号序列，Ni为第i次添加的白噪声序列。

3）对各个分量取平均，得到：

4）对n次叠加白噪声产生的分量对应相加，求其平均值，最终得到CEEMD分解结果：

其中，cj为CEEMD分解后的第j个本征模态分量。

CEEMD方法通过对原始信号叠加正负一对的高斯白噪声，增加了原始信号的频域范围，使得信号在各个尺度连续起来，在构建包络线时，极值点分布均匀，进而减少了模态混叠现象。正负一对的白噪声相互抵消，信号重构时对重构信号的影响可以忽略。

3 基于形态-CEEMD的地磁测量信号滤波方法

由于载体磁噪声和环境磁噪声种类复杂，为更好地保持目标信号的几何结构特征，采用广义形态滤波器对原始信号进行预处理，以滤除尖峰和正负脉冲。

3.1 结构元素的选取

结构元素包含三个特征量，即长度L，映射方式和幅值。不同的特征量对信号噪声抑制有着较大的影响。长度L往往大于噪声脉冲宽度，并小于信号长度，幅值一般不大于有用信号的1%，信号直线型和圆型结构元素能够有效滤除白噪声，三角型和正弦型对脉冲噪声的滤除效果较好。本文对比了不同长度下不同映射的滤波效果。

3.2 模态分量阈值选择方法

针对水下磁场环境与载体随机干扰磁场持续时间短、幅值大且伴随整个测量过程的特点，使用CEEMD方法对信号进行分解，对CEEMD分解后的固有模态分量（IMFs）逐一计算其自相关函数，通过计算自相关系数，分析各模态分量的自相关度确定阈值，选择相关度大于阈值的信号进行重构，并对舍弃信号进行均值滤波后与有用信号共同重构去干扰后的信号。降噪后的信号可以表示为

其中：a1.a2，…，am为自相关度大于阈值的m个模态分量的阶数，b1，b2，…，bk为噪声分量的阶数，Wj(t)为均值滤波后的噪声分量。

3.2.1 自相关函数与信号相关性

自相关函数体现出函数在不同时刻的相关程度。假设x(t)为随机信号，则其自相关函数为

其中，τ为信号偏移时间。在零偏移时，信号的自相关函数达到最大，理想高斯白噪声的归一化自相关函数在零偏移处值为1，其余偏移均为零。经过CEEMD分解后的噪声分量往往包含多个信号模态和非白噪声序列，通过观察不同信号的自相关函数可以看出，间断信号的自相关函数也是间断的，能量集中在零偏移附近。当信号是周期信号时，其自相关函数也为周期信号：

由于不同信号的能量存在较大偏差，自相关函数不能完全表现出信号的时间相关性。自相关函数在零偏移时刻达到最大值，为了比较多个信号的自相关性大小，需要对各个信号的自相关函数进行归一化处理，处理后的归一化自相关函数为

归一化自相关函数避免了信号能量对标准差的影响，能比较信号之间的自相关性大小。本文使用各个模态分量的归一化自相关函数标准差作为IMF的评价参数，可由如下公式计算得到：

经过阈值判断后，信号主导分量有IMFa1，IMFa2，...，IMFam，噪声主导分量有IMFb1，IMFb2，...，IMFbk，

3.2.2 阈值选择与信号重构

将各个IMF的归一化去零点自相关函数标准差作为判定信号主导分量和噪声主导分量的依据。设定阈值函数，当模态分量的归一化去零点自相关函数标准差大于等于设定阈值时，则为信号主导分量，反之为噪声主导分量。阈值选取规则如下：残余分量为r(t)。噪声主导分量经过均值平滑滤波后，得到Wb1，Wb2，...，Wbk。将信号主导分量与滤波后的噪声主导分量和残余分量进行信号重构，得到去噪后的信号：

3.3 形态滤波-CEEMD地磁信号降噪方法

CEEMD方法虽然能通过叠加白噪声的方式在一定程度上减轻模态混叠现象，但是当处理的信号含有突变和较大脉冲干扰时，仍会出现较大的模态混叠现象。本文提出了使用广义形态滤波对原始测量信号进行预处理，滤除尖端脉冲和突变信号，接着对预处理后的信号进行CEEMD方法分解，结合噪声的自相关特性，提出使用自相关特性的模态选择方法，分离出噪声分量以及有用信号分量，并对噪声分量进行二次滤波，与有用信号分量共同构建目标信号。具体实现过程为根据原始信号特点和目标信号要求，确定结构元素的长度L、映射关系g(m)和幅值k；对原始信号进行广义形态滤波，得到预处理后的信号；对预处理后的信号进行CEEMD分解，得到多个本征模态分量和残余分量；计算各个本征模态分量的自相关度，根据阈值选择方法确定有用模态分量和噪声模态分量；对噪声模态分量进行中值滤波处理；将滤波后的噪声分量与有用模态分量进行信号重构，得到去噪后的地磁信号。

4 磁力仪实测信号分析

实验中使用磁通门磁传感器进行磁信号采集，数据采样率为200Hz，量程在±96000nT，分辨率为1nT。实验分析了磁通门传感器在强磁干扰环境下测量数据和弱磁干扰的突然变向的单轴测量数据，分别如图1和图2所示。

图1 强磁干扰实测数据

图2 突然变向时测数据

从图中可以看出，强磁干扰下的测量信号出现严重的脉冲干扰，在突然转向的情况下，存在较大的脉冲干扰。

4.1 不同结构元素的滤波效果

结构元素作为广义形态滤波的关键因素，其参数的设定对滤波效果有着较大的影响。根据实验数据特点，首先选取直线型和圆型作为结构元素基本形状，其表达式分别如下：

其中，k，L，m∈R。对于离散信号，通过改变k，L和m的采样点，可以改变结构元素的幅值和点数。针对强磁干扰实测数据，其脉冲干扰持续宽度Lq多在100个采样点左右，而突然转向脉冲干扰持续宽度Lz在200个采样点左右。分别选取，50，100，Lz=25，100，150，k=1，0.5，0.25，得到以下实验结果。

图3 强磁干扰实测数据不同长度结构元素滤波效果对比

图4 突然变向实测数据不同长度结构元素滤波效果对比

从图中可以看出，结构元素的长度对滤波有着较大的影响，当结构元素长度与尖端脉冲长度比较接近时，能够有效滤除尖峰脉冲，但是当结构元素长度过小时，尖端脉冲不能完全被去除，当结构元素长度过长时，会出现过度平滑现象，使信号失去细节信息。结构元素的形状决定了尖端脉冲被消去的弥补信号的形态，由于本文使用了直线型，滤波后的信号多出现直角形状，在突然转向信号的滤波中表现良好，能有效滤除因转向产生的尖端脉冲。

4.2 不同模态选择方法的滤波效果

经过广义形态滤波后的信号仍然存在部分高频干扰，使用CEEMD对其进行模态分解，并提出基于自相关度的模态分量选择方法选择有用信号进行信号重构。以强磁干扰实测数据为例，使用L=50，k=0.5的形态滤波结果作为CEEMD输入，叠加噪声幅值为10nT，叠加10次正负白噪声。分解得到如图5和图6模态分量。

根据基于自相关度的模态选择方法，选择IMF7-IMF11和残余分量直接进行信号重构，对IMF1-IMF6进行均值滤波后进行信号重构，均值滤波窗口与结构元素长度一致，得到目标信号如图5所示。可以看出，基于自相关度的模态选择方法要比传统模态选择方法滤波效果好，有效区分了噪声分量和有用信号分量。图6表示仅使用CEEMD方法和形态-CEEMD方法的滤波效果，可以看出形态-CEEMD方法能够更有效地滤除尖端干扰脉冲。

图5 不同模态选择方法滤波结果对比

图6 形态-CEEMD方法与CEEMD方法对比效果

5 结语

本文主要介绍了广义形态滤波和CEEMD方法，分析了不同结构元素下广义形态滤波对不同信号的滤波效果，结合噪声信号的自相关函数特点，提出了基于自相关度的模态选择方法。使用实测地磁数据对方法进行了方法分析和验证，试验结果说明提出的形态-CEEMD方法能够有效滤除强干扰磁脉冲和转向情况下的磁脉冲干扰。结构元素的长度对滤波效果有着较大的影响，虽然事先分析了不同长度对滤波效果的影响，选择了合适的元素长度，但是当脉冲干扰更加复杂时，这种方法存在一定缺陷，需要考虑自适应的结构元素选取方法。