0 引言
弹道信号一般是描述飞行物体的运动状态(位置、速度、加速度、姿态、角速度等),带有噪声的、离散的测量数据或待估计的状态信息。原始的包含噪声的弹道信号不能直接用于精确导航、制导、控制、轨迹预测或数据分析。因此,需要设计合适的滤波算法对其进行降噪,得到更接近真实运动状态的估计值,使轨迹看起来更连续、合理。即使某些传感器数据暂时缺失或不可靠,依然可以通过模型和滤波器持续提供最佳的状态估计。
滤波算法主要的作用是去除不需要的频率成分、增强感兴趣的信号部分或者使得信号满足特定的频率响应。传统的滤波器如低通、高通、带通滤波器分别用于去除信号中的高频、低频和用于选择特定频率范围内的信号。自适应滤波器则能够根据信号的实时性质自动调整其滤波特性,且自适应滤波器对于处理非线性、时变的信号干扰问题非常有效。在实际应用中,通常无法得到信号和噪声的先验统计特性,而自适应滤波器则能够实时地调节当前时刻的滤波系数以适应随机信号的时变统计特性[1],可以实现一定意义下的最优滤波,因此更加适用于弹道输出信号的滤波处理。
最小均方(Least Mean Square,LMS)算法是一种基于维纳滤波理论,具有较为简单的结构和相对稳定的性能,同时其计算的复杂度也相对较低,更加易于硬件的实现,实际应用更为广泛。郜等[2]提出一种利用最小均方算法,快速估计出半球谐振子的特征参数,所提算法不仅可以对半球谐振子的性能进行量化评估,且进一步为半球谐振陀螺的误差分析与补偿提供理论依据。Rahaman等[3]在归一化阵列因子和均方误差的基础上测试了LMS算法在智能天线系统中的性能,测试结果显示LMS算法计算复杂度低,且在50次迭代后收敛。王等[4]基于LMS滤波算法优化数控加工误差控制,进一步提高了数控机床加工速度和加工精度,并有效解决了数控机床加工运动过程中存在的加工误差问题。张等[5]提出一种基于扩张状态观测器的永磁同步电机自适应逆控制方法,通过非线性滤波器以及归一化LMS算法实时估计转速,利用估计的权值对电机系统进行在线逆建模,研究表明控制方法具有更好的动态性能和抗扰能力。然而,传统的LMS自适应滤波器的稳态误差依然较大,同时该问题也是现阶段亟待突破的研究内容。
目前,学者们主要基于改变步长因子提出了一系列改进的LMS算法,从而优化该算法的收敛性能和稳态特性。Wang等[6]根据估计误差功率和系统噪声功率的平方根调整滤波器步长,提出了一种有效的非参数变步长归一化LMS算法,所提算法在各种信噪比情形下表现出优越的性能。SAEED等[7]提出了一种用于低信噪比应用的变步长增量LMS算法,该算法策略在低信噪比值下依然具有较高的性能。雷等[8]基于S函数的曲线特点,通过对函数的平移变换得到LMS算法步长因子的表达式,该算法在收敛速度、稳态误差及抗干扰能力方面的性能均有了很大的改善。刘和张等[9,10] 在对回波信号的滤波处理过程中,采用LMS自适应滤波算法进一步提升了滤波效果。吴等[11]提出一种改进的变步长LMS算法,该算法通过对传统的LMS算法中加入步长反馈因子并更新权值,使得所提算法兼顾了收敛速度和稳定性,且可以有效滤除噪声信号。张等[12]基于步长因子与误差信号的非线性函数提出了一种变步长的LMS自适应滤波算法,该算法兼顾了收敛速度和稳态误差,并将其成功应用于水位监测数据的滤波过程。叶等[13]基于双曲余弦函数,通过改变外部参数从而调节LMS滤波器步长因子,由此设计的自适应噪声对消系统,有效解决了引信信号处理系统中的噪声干扰问题。张等[14]在对无线电引信进行信号处理的过程中,提出了一种新的基于迭代步长、误差和动态约束迭代步长的数学模型,该模型算法具有很好的收敛性能和稳态特性,进一步提高了抗干扰能力。火等[15]采用类箕舌线函数代替Sigmoid函数作为步长迭代公式,提出了一种变步长归一化自适应滤波算法,极大的提高了算法稳定性。李、张和孔等[16-18]分别提出了基于反正切函数和双曲正切函数的变步长LMS自适应滤波算法,解决了收敛速度和稳态误差的矛盾。徐等[19]利用基于改进步长因子LMS算法的自适应滤波器对轴承-转子系统进行振动控制,实验结果表明该算法可以有效识别并滤除与转子同频的信号分量,从而实现储能飞轮主动磁轴承-转子系统的不平衡振动控制。冯等[20]提出一种基于变步长凸组合LMS自适应滤波动态血氧检测方法,对原始光电描记脉搏波信号经预处理后,将固步长与变步长LMS自适应算法进行组合,实验结果显示,该血氧检测方法对于高原旅游中高原反应的实时监测与提醒具有实际意义。张等[21]提出一种基于改进归一化最小均方(NLMS)的多雷达对消抗主瓣干扰方法,通过消散因子对步长实现自适应控制,实验结果表面该方法的收敛速度较快。徐等[22]构建了融合深度强化学习和LMS自适应滤波的算法模型,该算法不再使用传统的非线性函数调整步长,而是利用深度强化学习智能的控制步长因子,该算法在不同的初始条件下都能最快的调整至合理的步长值。
除了采用变步长的方法,也有学者们采用改变滤波器阶数来改善LMS自适应滤波算法的性能。孙等[23]提出了一种基于数字滤波器的自适应最优滤波器阶数的设计方法,实验结果表明,滤波器阶数是影响卫星导航接收机时域自适应抗干扰性能和计算复杂度的核心参数,所提方法通过自适应调节滤波器阶数,可以有效提升时域自适应滤波器的抗干扰性能。上述研究成果表明,滤波器的步长和阶数是优化LMS滤波算法的核心参数,因此有学者同时调节滤波器的步长和阶数,最大程度的优化LMS自适应滤波算法。冯等[24]通过粒子群优化算法自适应的调节LMS算法中的步长和阶数,并利用该算法设计了数字基带传输系统,使得该系统在噪声环境下具有更加优越的性能。但是,由于传统的粒子群优化算法在迭代寻优过程中容易陷入局部最优,因此获得的步长和阶数并不是全局最优值。本文提出一种新的LMS自适应滤波算法,通过采用改进的粒子群优化算法来自适应的调整LMS滤波器的阶数和步长,并且可以使得粒子群在迭代过程中获取到全局最优值,最大程度的优化LMS自适应滤波算法,并将所提算法应用于弹道信号的滤波过程。
1 改进的LMS自适应算法滤波
1.1 LMS自适应滤波原理
最小均方自适应滤波器包括FIR结构的维纳滤波器和LMS算法两部分[1]。首先,给定LMS自适应滤波器的M组滤波器的采样输入信号u(j-1),u(j-2),...,u(j-M),滤波器在时刻j的输出信号${\widehat{u}}_{j}$为输入信号的线性组合,也是当前时刻采样值uj的估计值:
${\widehat{u}}_{j}$=$\sum _{i=1}^{M}$ciu(j-i)=CTU
其中M为滤波器阶数,ci(i=1,2,...,M)为滤波器权系数,C为权系数ci的向量形式,U为输入信号的向量形式。时刻j的估计误差${f}_{{M}_{j}}$为
${f}_{{M}_{j}}$=Dj-${\widehat{u}}_{j}$(j=1,2,...,J)
其中,J为弹道信号总的时间点数,Dj表示理想输出信号。根据LMS算法可定义滤波器权系数的迭代公式:
ci=ci-1+2μ·fM·U
μ为LMS滤波器的迭代步长,误差向量fM=${\left[{f}_{{M}_{1}},{f}_{{M}_{2}},...,{f}_{{M}_{J}}\right]}^{T}$。LMS自适应滤波器的代价函数定义为
Fj=E[${f}_{{M}_{j}}^{2}$](j=1,2,...,J)
由公式(4)不难得出,${\widehat{u}}_{j}$的取值决定了代价函数的取值,而${\widehat{u}}_{j}$中滤波器阶数M和权系数ci的取值至关重要,且权系数ci的取值取决于滤波器的迭代步长μ。因此,LMS自适应滤波算法需要解决的问题是如何选取最优的权滤波器阶数M和迭代步长μ,使得估计误差的最小均方误差值最小。
1.2 IPSO-LMS自适应滤波算法
本小节给出了本文所提出的一种改进的LMS自适应弹道信号滤波算法的具体步骤。
步骤1:首先对粒子群智能优化算法进行初始化。给定LMS自适应滤波器的迭代步长、滤波器阶数的初始参数值,然后给定种群规模R,最大迭代次数G,惯性权重ω以及位置范围[Xmin,Xmax],学习因子l1=l2=1,l3=1.4。给定速度范围:
Vmax=Xmax/10=-Vmin
其中,初始速度设为0,初始位置为
${x}_{i}^{0}$=(Xmax-Xmin)·rand+Xmin
其中,rand表示0到1之间均匀分布的随机数。
步骤2:然后将初始产生的全部粒子分为三组,每组粒子采用不同的搜索机制。第一组粒子的搜索方向如下:
$\begin{array}{l}{v}_{i}^{k+1}=\omega {v}_{i}^{k}+{l}_{1}{r}_{1}(Pbes{t}_{i}^{k}-{x}_{i}^{k})\\ +{l}_{2}{r}_{2}(Gbes{t}^{k}-{x}_{i}^{k})+{l}_{3}{r}_{3}({U}^{k}-{x}_{i}^{k})\end{array}$
${x}_{i}^{k+1}$=${x}_{i}^{k}$+${v}_{i}^{k+1}$
其中,${v}_{i}^{k+1}$和${x}_{i}^{k+1}$分别表示第i个粒子在第k+1代的速度和位置,Pbes${t}_{i}^{k}$表示第i个粒子在第k代的个体最优位置,Gbes${t}_{}^{k}$表示整个粒子群体在第k代的群体最优位置。ω表示粒子更新的惯性权重值。r1,r2和r3为[0,1]之间的随机数。Uk的表达式为
Uk=$\sum _{t=1}^{\left[0.9\right.n]}$At/[0.9n]
其中,[·]表示取整,n表示高质量搜索集的元素数量。
第二组粒子的搜索方向如下:
${v}_{i}^{k+1}$=(${x}_{i}^{k}$-Gbestk)·rand·α
${x}_{i}^{k+1}$=${x}_{i}^{k}$+${v}_{i}^{k+1}$
其中,α=0.5表示速度缩放比例,用于控制粒子的飞行区域。第二组粒子在更新前先产生随机数R1,在R1>0.5的情况下根据公式(10)和公式(11)更新其位置,否则根据公式(7)和公式(8)更新其位置;第二组粒子的搜索机制是以当前全局最优粒子的反方向作为搜索速度,目的是进一步增加粒子的多样性。
第三组粒子的搜索方向如下:
${x}_{i}^{k+1}$=$\sum _{t=1}^{\left[0.5n\right]}$Bt/[0.5n]
其中,Bt为集合B中的元素,B⊆S,B中元素的数量为[0.5n]。对于第三组粒子,在更新前先产生随机数R2,在R2>0.5的情况下根据式(12)更新其位置,否则根据公式(7)和公式(8)更新其位置。第二和第三组粒子是对种群的变异,子种群均取5,变异概率为0.5。
步骤3:根据上述搜索机制分别更新LMS自适应滤波器的迭代步长和滤波器阶数,并根据LMS自适应滤波器的输出信号和理想输出信号定义改进的粒子群智能优化算法的适应度函数。
本算法中,每组粒子包含两维变量,因此i=1,2。第k代粒子可表示为${x}_{1}^{k}$和${x}_{2}^{k}$,其中${x}_{1}^{k}$表示LMS滤波器阶数,${x}_{2}^{k}$表示LMS滤波器的迭代步长。值得注意的是,滤波器阶数取正整数,但粒子群算法中寻找的${x}_{1}^{k}$的最优粒子有可能不是正整数。因此,取${x}_{1}^{k}$=2·round(${x}_{1}^{k}$/2),“round”表示四舍五入,与matlab中的“round函数”含义一致。在时刻j,LMS滤波器的输出信号Yj为
Yj=$\sum _{i=1}^{{x}_{1}^{k}}$ciu(j-i)=CTU
滤波器权系数的迭代公式为
ci=ci-1+2·${x}_{2}^{k}$·error·U
误差向量error=${\left[{e}_{1},{e}_{2},...,{e}_{J}\right]}^{T}$,ej表示理想输出信号与滤波输出信号之间的误差:
ej=Dj-Yj(j=1,2,...,J)
其中,J为弹道信号总的时间点数,Dj表示理想输出信号。
第k代粒子的适应度函数为
minX∈Ω$\sum _{j=1}^{J}\Vert {Y}_{j}({x}_{1}^{k},{x}_{2}^{k})-{D}_{j}{\Vert }^{2}$
其中,X={${x}_{1}^{k}$,${x}_{2}^{k}$},k=1,...,G,G为最大迭代次数,Ω为所有粒子构成的粒子群空间。
步骤4:根据所有粒子的适应度函数值更新个体极值和群体极值。设置一个小量h=0.01,若某一粒子的适应度函数值小于h,就把该粒子加入搜索集,否则就舍弃该粒子。若搜索集中的粒子数量大于500,就仅取500个粒子,且这500个粒子的适应度值最小,否则就保留搜索集中全部的粒子。判断迭代次数k是否小于最大迭代次数G,若小于则令k=k+1,并回到步骤2;否则输出全局最优解Gbest。
根据上述步骤获得LMS算法的最优迭代步长和最优滤波器的阶数,在此基础上根据公式(13)-(15)计算得到LMS滤波器对弹道信号滤波后的输出信号Y*。图1为基于改进的粒子群优化算法和最小均方自适应算法(IPSO-LMS)的信号滤波方法的流程图。
2 算例分析
2.1 IPSO-LMS与PSO-LMS算法的滤波效果对比
首先由导弹测控系统测得含有噪声干扰的弹道信号,然后通过理论模型推导得到的数据作为理想输出信号。本节选取理想输出信号为参考信号,将已有的标准粒子群优化算法与LMS算法相结合(简称PSO-LMS)滤波器的输出信号作为对比信号,来验证本文所提方法IPSO-LMS滤波器的性能。图2中红色线条为某弹道的理想输出信号,蓝色线条为PSO-LMS滤波器的输出信号,黑色线条为本文所提方法IPSO-LMS滤波器的输出信号。图3中蓝色信号为PSO-LMS滤波器的输出信号与理想输出信号的偏差,红色线条为IPSO-LMS滤波器的输出信号与理想输出信号的偏差。由图2和图3可以看出,IPSO-LMS滤波算法得到的输出信号更加平滑,稳态响应的收敛速度更快,且稳态误差几乎为零。
2.2 不同SNR下IPSO-LMS算法的滤波效果
2.2.1 滤波效果时域对比
图6 弹道滤波信号的时域对比(SNR=0 dB)Fig.6 Comparison of ballistic filtered signals in time domain (SNR=0 dB) |
图7 弹道滤波信号的时域对比(SNR=5 dB)Fig.7 Comparison of ballistic filtered signals in time domain (SNR=5 dB) |
图8 弹道滤波信号的时域对比(SNR=10 dB)Fig.8 Comparison of ballistic filtered signals in time domain (SNR=10 dB) |
2.2.2 滤波效果频域对比
图10到图13分别给出了原始信号(红色线条)、加噪信号(蓝色线条)以及IPSO-LMS算法的滤波信号(黑色线条)的频域对比图,其中信噪比(SNR)分别取0,5,10和20 dB。由图10可以看出,在SNR=0 dB时,加噪曲线(蓝色线条)在全频段有明显噪声波动,滤波曲线(黑色线条)显示出噪声能量大幅降低,说明IPSO-LMS算法能有效压制噪声。由图10到图13可以看出,SNR越小(噪声强),滤波后的噪声能量削弱越显著,噪声抑制效果越好。同时,随着SNR越大(噪声弱),滤波曲线与原始曲线的贴合程度越高。当SNR=20 dB时,由图13可以看出滤波曲线和原始曲线在低频段几乎重合,说明滤波信号对原始信号的能量保留较好。
图10 弹道滤波信号的频域对比(SNR=0 dB)Fig.10 Comparison of ballistic filtered signals in frequency domain (SNR=0 dB) |
图11 弹道滤波信号的频域对比(SNR=5 dB)Fig.11 Comparison of ballistic filtered signals in frequency domain (SNR=5 dB) |
图12 弹道滤波信号的频域对比(SNR=10 dB)Fig.12 Comparison of ballistic filtered signals in frequency domain (SNR=10 dB) |
2.2.3 滤波性能指标分析
3 结论
针对弹道信号中存在的噪声干扰问题,尤其是没有信号的先验统计特性,因此不能使用某一固定参数的滤波器来处理该噪声干扰问题。本文基于改进的粒子群优化算法,提出一种改进的最小均方自适应的信号滤波方法。由于迭代步长、滤波器的阶数决定了LMS自适应滤波器的性能,因此合理选择LMS算法的迭代步长和滤波器的阶数至关重要。传统的LMS自适应滤波器一般采取遍历的形式选择迭代步长、滤波器的阶数,这是十分耗时的,而且不能保证选取到最优性能的参数。此外,传统的标准粒子群优化算法中,所有的粒子均采用相同的更新策略,导致种群多样性缺失,并且标准粒子群优化算法还特别容易陷入局部最优。本文进一步改进粒子群优化算法,通过智能的筛选出迭代步长、滤波器阶数的最优参数值,从而使得LMS自适应滤波器具有最优的性能。智能的筛选出LMS滤波器的阶数和迭代步长的最优参数值,进一步提升了LMS自适应滤波器的性能。数值结果表明,与PSO-LMS算法相比,本文所提的IPSO-LMS算法输出信号更加平滑,前2秒内瞬态响应振荡幅值大幅度减小,稳态响应的收敛速度更快。此外,在不同输入SNR下,本文所提的IPSO-LMS滤波算法对于低SNR信号(噪声强)的降噪效果更加突出,即在一定范围内大幅提升SNR可以一定程度减小均方根误差。对于高SNR信号(噪声弱),虽然均方根误差进一步缩小,但SNR优化空间有限,因此本文所提的IPSO-LMS滤波算法在低SNR下滤波效果优势明显。