连续变速颤振试验的自适应粒子滤波算法

谭博

(航空工业第一飞机设计研究院 综合航电系统设计研究所， 西安 710089)

0 引 言

颤振试验数据处理的基本目的是分析实测亚临界响应信号，完成模态参数估计(Modal Parameter Estimation，简称MPE)和颤振边界预测(Flutter Boundary Prediction，简称FBP)。目前常用的风洞颤振试验方法为台阶法，即在一个台阶时间内保持外界条件基本稳定并采集振动信号，因信号模态参数在台阶时间内变化很小，可认为采集的信号是近似平稳的。通过连续采集多个台阶的振动信号并提取其模态参数，以阶梯的方式预估信号模态参数变化的趋势，并达到预测颤振边界的目标。此方法易于操作，但因其与颤振发生的实际原理不符，导致预测结果存在较大的偏差。为了解决台阶试验方法的问题，一种连续变速颤振试验(Flutter Test with Variable Progression Speed，简称FTVPS)方法被提出，该方法是马赫数、速度、高度等飞行参数个别或全部随时间不断变化的一类颤振试验。与传统的台阶法相比，其试验周期短且更符合实际状态。但是，试验采集的信号为非平稳信号，在亚临界状态下，各项特征变化剧烈，为试验数据的处理和分析带来了新的困难。

非平稳信号广泛存在于各工程应用领域[1]，且各领域的非平稳信号之间存在很大的差异性，这使得非平稳信号的分析处理方法通常具有很强的专业领域限制和问题针对性。例如，用于分析切削振动信号的经验模态分解方法[2]、时域和频域分开的非平稳随机振动分析方法[3]以及基于最大谱的非平稳随机振动数据分析处理方法[4]等就是针对特定工程领域问题提出的方法，这些方法一般不具有普适性，且多用于处理具有随机激励的振动信号，而颤振属于自激振动，因此上述方法难以在连续变速颤振试验数据处理中应用。

连续变速颤振试验数据处理一般使用时变参数建模方法，主要有自适应滤波、基函数展开和粒子滤波三种。其中，自适应滤波算法和基函数展开算法在处理信号时，容易受到噪声的影响[5]，因此连续变速颤振试验数据处理中一般使用粒子滤波算法。粒子滤波算法是20世纪50年代末提出的一种基于贝叶斯采样估计的顺序重要采样滤波方法[6]，此后，N.Gordon等[7-8]于1993年又提出了新的基于顺序重要采样的Bootstrap非线性滤波算法，该方法的目标是完成目标状态的在线跟踪。目前，在各工程领域内的非平稳信号处理中有广泛应用，例如，王宏健等[9]提出一种改进差分粒子滤波算法；陈志敏等[10]提出一种自适应粒子群优化的目标跟踪算法；叶华等[11]提出一种基于稀疏表示的粒子滤波算法；林晓梦等[12]提出改进粒子滤波重采样算法；谌剑等[13]进行了粒子滤波算法的权值优化组合研究；郑华等[14]提出一种用于多传感器组合测量的粒子滤波算法。以上粒子滤波算法都针对应用时面临的问题进行了改进。

本文针对连续变速颤振信号的特点并结合工程实际，选用合适的非平稳度的度量方法，并将非平稳度引入粒子滤波算法中，并通过仿真实验验证该方法分析非平稳程度变化信号的有效性。

1 非平稳度

平稳信号一般指信号特性时不变的信号，与之对应的非平稳信号则是指信号特性会随时间变化的信号。非平稳信号种类多样且都具有各自的特点，因此很难找到一种可广泛适用的方法度量非平稳信号特性变化剧烈程度(Non-stationary Degree，简称ND)。目前存在的度量方法有很多，如周期非平稳度(Degree of Cyclostationary,简称DCS)，基于Hilbert时频谱的ND，基于相关积分值的ND等。上述三种度量方法中，DCS一般针对周期非平稳信号，对非周期信号不具有适用性；第二种计算复杂，且其与平稳信号的关系不够明确；第三种针对具体的工程问题，均不适用于处理本文研究的问题。

连续变速颤振试验中测量信号一般为非周期信号，且非平稳程度会随时间发生变化，在亚临界状态下变化速度明显增大。因此，本文选用文献[15]中给出的非平稳度计算方法，该方法得到的ND值恒大于0，其中，平稳信号的ND值小于1，非平稳信号的ND值一般大于1且与其幅值或频率变化速率呈近似线性关系。

2 自适应粒子滤波算法

对于连续变速颤振信号这一非线性信号，数学推导往往无法进行，一般使用粒子滤波算法结合蒙特卡洛方法进行估计。具体方法为：以一组随机或已知的某种分布样本来描述估计量的概率分布，再根据得到的测量值，通过重要性函数对各样本点的权值进行调整，以该带权值的样本序列来逼近真实的后验概率分布，从而序贯更新状态，其核心步骤是通过重要性函数进行的采样和调整过程，即重采样过程。

连续变速颤振试验信号的一个主要特点是信号的幅频特性随时间非线性变化，在初始阶段变化缓慢，在临界状态会剧烈变化直至颤振发生结构破坏。这使得一般的粒子滤波算法在采集信号的幅频特性变化速度增快后，无法进行有效的跟踪，导致最终处理结果和实际结果相差过大，无法准确预估颤振边界。因此，需要针对其非平稳程度变化的特点优化调整重采样过程，将非平稳度作为重要性函数的一个变量进行综合考虑。

考虑到非平稳度增加后，需要在更广泛的状态空间内寻找权值更大的粒子，本文的自适应滤波算法在一般的重采样过程后，引入粒子调整步骤，依据非平稳度的变化，增加或减少粒子数量及步长，达到扩大或缩减粒子群在状态空间内分布的目的。同时，为了避免非平稳度过大引起粒子数剧增，导致无法计算，或非平稳度过小引起粒子数剧减，导致粒子数不足的情况，对步长和粒子数调整幅度进行限制。最终可得到如下自适应粒子滤波算法步骤。

初始化：t=0

当t=1,2,…时，重复进行以下步骤：

(1) 重要性采样

(2) 重采样

(3) 调整粒子

计算当前ND值，并依据ND值调整粒子群：

其中，Pt′为调整粒子集；μ为调整的步长；ND为信号当前的非平稳度；μmax为可调整的最大步长；μmin为可调整的最小步长。

最终期望：

可以看出：当信号ND值增大时，本文方法将对粒子集Pt进行扩充，扩充后的新粒子集粒子数量更多，在状态空间中的分布更加分散，可在更大范围内寻找更精确的状态；当信号ND值减小时，本文方法将缩减粒子集以提高粒子的聚集程度，减少运算量并提高运算速度。

3 仿真试验

生成两组变ND信号，其ND值由初始值起，于第2、第4和第6 s时分别增加为初始值的2倍、3倍和4倍。信号时长8 s，采样频率为128 Hz。生成信号在无噪声环境下的时间历程及经由STFT求得的时频分布分别如图1～图2所示。向两个信号中依次混入信噪比不同的噪声后得到两组仿真信号。

图1 变频变ND信号时频分布

图2 变幅变ND信号时频分布

分别使用一般粒子滤波算法和本文的自适应粒子滤波算法对两组生成信号进行处理，所得的跟踪相对误差结果如图3～图6所示，由上至下依次为无噪声和信噪比分别是10、6、3、0 dB。

图3 一般粒子滤波算法(PF)变频信号跟踪结果

图4 一般粒子滤波算法(PF)变幅信号跟踪结果

图5 自适应粒子滤波算法(APF)变频信号跟踪结果

图6 自适应粒子滤波算法(APF)变幅信号跟踪结果

对比图3和图5，可以看出：对于变频变ND信号，无噪声条件下，本文方法跟踪精度高于一般粒子滤波算法，其误差约为一般粒子滤波算法的50%；有噪声时，本文方法跟踪精度与一般粒子滤波算法相近。

对比图4和图6，可以看出：对于变幅变ND信号，本文方法的跟踪精度明显高于一般粒子滤波算法，其误差约为一般粒子滤波算法的10%～25%。

由实验结果可以得到如下结论：

(1) 对于变频信号，自适应滤波算法在无噪声条件的跟踪精度高于一般粒子滤波算法；

(2) 变频变ND信号在ND值发生变化的瞬间，自适应滤波算法的误差存在跳变情况；

(3) 对低信噪比的变频信号，自适应滤波算法与一般粒子滤波算法跟踪精度相近；

(4) 自适应粒子滤波算法的跟踪精度随着ND值的增大而降低，随着信噪比的降低而降低；

(5) 自适应粒子滤波算法对变幅信号的跟踪精度优于一般粒子滤波算法。

4 结 论

本文提出了一种依据信号非平稳度来调整重采样函数的自适应粒子滤波算法，对变幅变非平稳度信号的跟踪精度优于一般粒子滤波算法，且在无噪声条件下，对变频变非平稳度信号的跟踪也具有更高的精度，为连续变速试验数据分析打下了基础，也为其他领域的非平稳信号分析处理方法提供了一种解决思路。