奇异谱分解在超声速无人机声振试验数据处理中的应用

第一作者刘鎏男，博士生，1986年生

通信作者闫云聚男，教授，1954年生

奇异谱分解在超声速无人机声振试验数据处理中的应用

刘鎏，闫云聚，李鹏博

(西北工业大学力学与土木建筑学院，西安710072)

摘要：将相空间重构和奇异谱分解相结合对受强气动噪声影响的超声速飞行器测试数据进行滤波,以实现对于试验参数的精确识别。首先通过数值仿真论证该方法的可行性,然后针对某型超声速无人机的声振试验对试验采集数据进行相空间重构,并对重构后的轨迹矩阵进行奇异值分解,得到反映真实信号信息的信号子空间和反映噪声信息的噪声子空间。通过定义奇异值差分谱这一指标来判定真实信号信息子空间维数,并针对现有最大差分谱理论缺陷提出了优选差分谱峰值理论,利用奇异谱分解的逆过程对真实信号进行重构。重构结果表明,该方法适用于超声速飞行下的飞行器声振试验数据处理,为超声速飞行器飞行状态的精确描述提供了良好的思路。

关键词：相空间重构；奇异值分解；超声速；声振试验；信号滤波

收稿日期：2013-12-11修改稿收到日期：2014-01-11

中图分类号:TB53文献标志码：A

Singular value spectral decomposition and its application in acoustic vibration test data processing of a supersonic aircraft

LIULiu,YANYun-ju,LIPeng-bo(College of Civil Engineering and Mechanics, Northwestern Polytechnical University, Xi’an 710072, China)

Abstract:In order to realize precise identification of acoustic vibration test data of a supersonic aircraft, a new filtering method combining phase space reconstruction and singular spectral decomposition was proposed. Firstly, the feasibility of this method was demonstrated through numerical simulation. Secondly, in order to separate the signal subspace and the noise subspace, the phase space reconstruction of the test data was conducted, and the attractor track matrix was also decomposed with singular value decomposition (SVD). Finally, aiming at shortages of the maximum difference spectrum theory, the concept of optimizing difference spectrum theory was presented, and the signal reconstruction was proposed on the basis of the peak position of the optimizing difference spectrum. Reconstruction results showed that the proposed method is suitable for processing the acoustic vibration test data of a supersonic aircraft, the result provided a good foundation for the precise description of a supersonic aircraft’s flying state.

Key words:phase space reconstruction; singular value decomposition; supersonic; acoustic vibration test; signal filtering

超声速飞行器在高速飞行时,所形成的复杂气动流场会产生强大的气动噪声,并且有可能与飞行器结构体形成严重的耦合振动效应。飞行器结构壳体外面的强大气动噪声和其自身的强烈振动都会传递到飞行器仪器设备舱内,使内部的仪器设备处于强烈的振动与噪声环境中,对设备造成的可能危害包括：导线擦伤、紧固件松动、间断电触点、带电元件的接触和短路、密封失效、构件疲劳、光学上的失调、裂纹和断裂、微电子元件线路故障、导线的磨损、印刷电路板破裂等各种故障,将严重影响飞行器飞行的可靠性与安全性。为了预防以上可能的危害并进一步准确预示飞行器仪器设备的振动噪声环境,需要进行声振试验从而对噪声环境下的飞行器振动响应进行预报。试验中振动响应数据的获取不可避免地存在各种干扰,造成诸如测试数据模糊,信噪比低,信号突变等问题,直接影响到后续的数据处理和试验理论论证对比分析。因此,如何选取适当的信号滤波方法对试验数据进行预处理就显得尤为重要。

近几年来奇异谱分解技术(Singular Value Decomposition, SVD)被广泛运用。张友民等[1-2]和宫晓琳等[3]将基于SVD的固定区间平滑方法应用于机载合成孔径雷达运动补偿用位置姿态系统的后处理,获得了位置、速度和姿态的最优估计;Shim等[4]和贾鹏等[5]将基于SVD的冗余惯导系统故障诊断方法用于双陀螺故障检测;Shnayderman等[6]和朱卫纲等[7]提出了基于SVD的遥感图像融合性能评价,较好地反映了多类遥感图像融合的质量。何田等[8-10]将SVD运用于噪声背景下突变信息的检测,实现了对于某型航空发动机碰摩故障实测数据的处理。

SVD是一种非线性的滤波方法,它克服了傅里叶变换只能处理平稳信号的局限性，不需要像小波分析那样耗时耗力的凭借经验选择小波函数,也没有Hilbert-Huang变换所存在的边界处理以及模态混叠等问题，与上述方法相比,它具有滤波算法简单,滤波能力强以及频域范围普适性等优点。过往的一些关于SVD方法的研究往往集中于对一些中低频信号进行研究,信号本身较为简单且所处环境亦较为单一,少有关于SVD方法运用于复杂环境下振动测试数据处理的相关研究。鉴于此,本文在以往研究的基础上,将SVD方法和相空间重构相结合对超声速飞行器声振试验数据进行处理。

1相空间重构与奇异值分解

对于一系列离散时间序列{hi,i=1,2,…,N},其中N为采样点编号,利用Takens定理[11]可以对其进行相空间重构,得到如下L×K阶Hankel矩阵：

(1)

其中：K=N-L+1,延时值为1,H也被称为轨迹矩阵。

轨迹矩阵H经SVD分解后可得如下形式：

H=UΣVT

(2)

式中：U∈RL×L，Σ是一个L×K的对角阵,Σ=diag(σ1,σ2,…,σn)，σ1,σ2,…,σn为奇异值且按降序排列,V∈RK×K。

式(2)两边同时右乘V,有

HV=UΣ

(3)

当离散时间序列{hi,i=1,2,…,N}含有噪声时,轨迹矩阵H作SVD分解后可表示为：

UsΣsVTs+UωΣωVTω

(4)

由式(4)可以看出,轨迹矩阵经奇异值分解以后,得到了表征信号信息的信号子空间H′和表征噪声信息的噪声子空间N,这两个空间是互不相关的,SVD分解的目的就是将这两个互不相关的子空间分离以达到信号滤波的作用。

2相空间重构维数和最小嵌入维数确定

2.1相空间重构维数确定

前文已经提到,只有在由一维观测序列及其适当延时值所构成的维数合适的相空间中,系统演化的动力学行为才能由重构空间中点的演化轨迹无奇异的表达出来。这里存在一个相空间重构维数的确定问题,即轨迹矩阵H的维数L×K的确定。徐锋等[13]通过反复试验的方法近似给出了相应确定原则而并未从理论上给予相应的证明，本文通过理论分析验证该原则的正确性。

(5)

定理2.1 对于式(1)所示Hankel矩阵H，矩阵的迹可以表示为

(6)

证明依据式(2) Hankel矩阵H的定义可以得到

(7)

上式可变换为

(8)

(9)

(10)

(11)

上式表明：

(12)

L=[(N+1)/2]

(13)

K=N-L+1

(14)

其中符号[]表示取整。

2.2最小嵌入维数确定

对于含有噪声的离散时间序列,经由奇异值分解后得到反映真实信号的信号子空间和反映噪声信号的噪声子空间,这里就存在两个子空间维数的确定问题。反映到实际运算过程中就是有效奇异值阶数选择问题。原则上,奇异值分解的滤波效果主要取决于有效奇异值的选择。当所选的奇异值数目较少,滤波阶次较低时,会导致滤波信号包含的信息不完整,甚至发生波形畸变现象,难以对原信号的有效信息特征做出准确反映;而当所选奇异值数目较多,滤波阶次较高时,在滤波后的信号中仍然保留了一部分噪声信息,无法达到信号滤波的目的。事实上,由式(4)可知,轨迹矩阵H作SVD分解后,所得到的奇异值矩阵具有如下形式：

(15)

其中Σs是反映真实信号信息的奇异值部分,而Σω则代表噪声信息的奇异值部分,即轨迹矩阵H作SVD分解后会分为反映真实信号信息的Σs部分以及代表噪声信息的Σω部分。由于真实信号信息和噪声信息之间存在差异性,反映到奇异值矩阵上在交接部分奇异值大小必然会存在突变,那么在Σω和Σs过渡处的奇异值在数值上会明显有一个突变点,也就是说突变点前的奇异值所对应的分量为真实信号信息,突变点后的奇异值所对应的分量为噪声信息。这里通过定义奇异值差分谱来衡量相邻奇异值的变化情况,即

(16)

其中：σi和σi+1表示相邻两奇异值数值大小,赵学智等[14-15]在对机械信号处理时采用了最大差分谱理论,即差分谱最大峰值处即为有效奇异值选择点。数值仿真表明,在外部干扰和实际信号差异性不大的情况下最大差分谱理论有较好的效果,但对于实际超声速飞行环境下的飞行器而言,其外部环境更为复杂,最大差分谱理论并不适用,这一点在后文中可见。为此,本文提出了优选差分谱峰值方法,即逐次找出奇异值差分谱上各个峰值点并分别进行信号重构,通过对比分析来寻找最优的重构曲线。

3滤波后信号的重构

最小嵌入维数即有效奇异值个数选定以后,下一步就是依据有效奇异值的个数来对滤波信号进行重构，当选定的有效奇异值个数为k时有

(17)

需要注意的是,H′并不具有Hankel矩阵的形式,它只是当选定矩阵秩为k时的一种最佳逼近,为了得到真实的信号序列sk,需要对矩阵H′中的反对角元素取平均,即

(18)

其中：p=max(1,i-K+1)，q=min(L,i)。

4数值仿真

为说明奇异谱分解的滤波效果以及验证重构相空间维数和最小嵌入维数的选择原则,构造信号

x(t)=sin(48πt)+sin(100πt)

(19)

在此基础上叠加符合正态分布N(0,1)的高斯白噪声n(t),噪声强度D=2。叠加噪声后信号可以表示为

y(t)=x(t)+n(t)=

sin(48πt)+sin(100πt)+n(t)

(20)

采样频率fs=1 000 Hz,采样点数1 024,原始信号和叠加噪声后信号见图1所示。

图1　原始信号与叠加噪声后信号 Fig.1 Original signal and noisy signal

图1(a)所示是加噪前的信号时域信息,图1(b)是加噪后的信号时域信息,可以看到加噪后信号波形无论在幅值还是形状上都变得杂乱无章，信号信息彻底淹没在噪声环境下。图1(c)和图1(d)分别表示加噪前、后的频域信息,同样可以看到加噪后频域部分会出现其他峰值。下面利用奇异谱分解技术对加噪信号进行处理,采样点数1 024,即N=1 024,由式(14)和式(15)有L=[(N+1)/2]=512，K=N-L+1=513,即相空间重构后轨迹矩阵H的维数为513×512,以此轨迹矩阵进行奇异值分解,得到的前21个奇异值见表1,相应的奇异值差分谱值如图2所示,从图2中可以看出奇异值峰值主要集中于前面部分,在i处于2,4,6,8,16,20和48时奇异值大小发生突变，且在i=4处出现最大峰值,δmax=122。利用式(17)和式(18)取前4个奇异值进行信号重构,重构结果如图3所示,可以看到时域方面两者基本吻合,重构效果明显,经过奇异值分解重构后信号的两个主要频率成分f1=24 Hz及f2=50 Hz得到很好地体现。可见奇异谱分解技术在信号滤波方面是切实可行的,同时奇异值差分谱值这一指标也能较准确反映有效奇异值的个数。

表1　前21个奇异值数值情况

图2　奇异值差分谱曲线 Fig.2 Singular value difference spectrum curve

图3　前4个奇异值分量重构结果 Fig.3 The first four singular value components reconstruction results

5奇异谱分解在超声速飞行器声振试验数据处理中的应用

通过声振试验对某超声速无人机的飞行状态进行描述，试验模型用弹性绳吊装于高声强混响室，如图4(a)所示，试验总声压级共五个状态，分别是总声压级140 dB、145 dB、150 dB、155 dB、160 dB，每个状态持续30 s，在飞行器表面，安装16个单轴向加速度传感器，加速度测点布局如图4(b)所示,加速度采样频率为64 kHz。不失一般性,本文取总声压级160 dB的测试数据进行分析,由于每个通道采集数据量较大,这里取其中的10 240个数据点进行分析,得到试验后测得的总声压级为160 dB时测点2的时频域信息如图5所示。

图5　加速度测点2#时、频域信息 Fig.5 Time domain and frequency domain of the acceleration data about sensor 2

下面利用奇异谱分解技术对实验数据进行处理。采样点数10 240,N=10 240,由式(14)和式(15)有L=[(N+1)/2]=4 096,K=N-L+1=4 097,即相空间重构后轨迹矩阵H的维数为4097×4096,以此轨迹矩阵进行奇异值分解,得到的前27个奇异值见表2,对应的奇异值差分谱曲线如图6所示。

表2　前27个奇异值数值情况

图6　奇异值差分谱曲线 Fig.6 Singular value difference spectrum curve

可以看到奇异值序列在i处于2,4,6,12,14,20,22时产生突变,且在i=2处出现最大峰值,奇异值差分谱最大值δmax=29 195,取前2个奇异值进行信号重构,重构信息如图7所示。

图7　前2个奇异值分量重构结果 Fig.7 The first 2 singular value component reconstruction results

图7表示的是重构后时、频域信息,可以看到在频域图中仅在218.01Hz处存在一个频率峰值,而噪声致振环境下结构的响应具有多频的特性[16]。这说明仅按照最大差分谱理论难以对超声速环境下飞行器声振响应进行准确描述。奇异谱分解的核心思想在于真实信号和噪声信号的分离,由于信号和噪声的差异性此种分离在临界处必然伴随着突变,奇异值突变处正是这种差异性的潜在体现。尽管高速飞行器所处振动环境复杂，声振信号存在多频的特性，但信号和噪声之间的差异性仍然存在。因此,这里同时考虑其它奇异值突变处的重构情况,重构结果见图8(a)～8(f)。

图8　重构结果 Fig.8 Reconstruction results

对比图1及图8(1)-8(6)可以看到,当取22个奇异值进行重构时效果最为显著,响应信号的主要频率成分见表3所示。

表3　声振信号主要频率成分

为了定量比较最大差分谱峰值理论和优选差分谱峰值理论，这里引入信噪比指标

(21)

图9　不同奇异值重构信噪比值 Fig.9 Signal to noise ratio under different singular value decomposition

可以看到按照最大差分谱理论σ=2时重构后信噪比SNR为2.30 dB,在所有峰值处奇异值重构结果中最小,当按照最优差分谱峰值理论选定有效奇异值σ=22时重构效果最为显著,信噪比SNR达到8.49 dB。

6结论

本文利用SVD和相空间重构实现了对于高速飞行器声振响应的状态描述,从理论上分析了相空间重构维数及有效奇异值阶数等指标并给出了相应的表达式及判定方法,为实际应用的参数选择提供了依据,并通过仿真分析和试验结果论证验证了该方法的可行性，结论如下：

(1)针对以往依据经验选择相空间重构维数的方法从理论上进行了分析论证，分析结果表明Hankel矩阵维数接近方阵时效果最佳；

(2)奇异谱分解以信号和噪声的差异性作为分界条件,通过寻找最优奇异谱峰值点可以确定有效奇异值个数；

(3)最大峰值理论并不适用于高速飞行器声振响应信号的信号重构，依据本文提出的优选差分谱峰值理论进行重构维数的选择效果较好,采用前者重构后所得信噪比仅为2.30 dB，而采用后者重构后信噪比能达到8.49 dB。

参考文献

[1]Zhang Y M, Dai G Z, Zhang H C, et al. A SVD-based extended Kalman filter and applications to aircraft flight state and parameter estimation[C]. Proceedings of 1994 American Control Conference, Maryland: Baltimore, 1994: 1809-1813.

[2]张友民, 焦凌云, 陈洪亮,等. 基于奇异值分解的固定区间平滑新方法[J]. 控制理论与应用, 1997, 14(4): 579-583.

ZHANG You-min, JIAO Lin-yun,CHEN Hong-liang, et al. A new fixed-interval smoothing algorithm based on singular value decomposition[J]. Control Theory and Applications, 1997,14(4): 579-583.

[3]宫晓琳, 房建成. 基于SVD的R-T-S最优平滑在机载SAR运动补偿POS系统中的应用[J]. 航空学报, 2009,30(2): 311-318.

GONG Xiao-lin, FANG Jian-cheng. Application of SVD-based R-T-S optimal smoothing algorithm to POS for airborne SAR motion compensation[J]. Acta Aeronoutica et Astronautica Sinica, 2009, 30(2):311-318.

[4]Shim Duk-sun, Yang Cheol-kwan. Geometric FDI based on SVD for redundant inertial sensor systems[J]. 2004 5th Asian Control Conference, 2004, 6:1094-1100.

[5]贾鹏, 张洪钺. 基于奇异值分解的冗余惯导系统故障诊断[J]. 宇航学报, 2006, 27(5): 1076-1080.

JIA Peng, ZHANG Hong-yue. Fault detection for redundant inertial navigation system based on SVD method[J]. Journal of Astronautics, 2006, 27(5):1076-1080.

[6]Shnayderman A, Gusev A, Eskicioglu A M. A SVD based grayscale image quality measure for local and global assessment [J]. IEEE Transactions on Image Process, 2005, 14 (2): 422-429.

[7]朱卫纲, 周萌清, 徐华平,等. 基于奇异值分解的遥感图像融合性能评估[J]. 北京航空航天大学学报，2008, 34(12): 1448-1451.

ZHU Wei-gnag, ZHOU Meng-qin, XU Hua-ping,et al. Remote sensing image fusion assessment based onSVD[J]. Journal of Beijing University of Aeronautics and Astronautic, 2008，34(12): 1448-1451.

[8]何田, 刘献栋, 李其汉. 噪声背景下检测突变信息的奇异值分解技术[J]. 振动工程学报, 2006, 19(3): 399-403.

HE Tian, LIU Xian-dong, LI Qi-han. An improved method of detecting abrupt information based on singularity value decomposition in noisebackground[J]. Journal of Vibration Engineering, 2006, 19(3): 399-403.

[9]何田, 刘献栋, 李其汉. 一种改进的航空发动机转静件碰摩故障诊断方法[J]. 航空动力学报, 2008, 23(6): 1093-1097.

HE Tian, LIU Xian-dong, LI Qi-han. Improved fault diagnosis method for aero engine rotor-sator rub[J]. Journal of Aerospace Power, 2008, 23(6): 1093-1097.

[10]何田, 刘献栋, 陈亚农,等. 基于奇异值差的转静件碰摩故障诊断方法研究[J]. 航空动力学报, 2009, 24(10): 2296-2301.

HE Tian, LIU Xian-dong, CHEN Ya-nong,et al. Method for diagnosing rub fault of rotor-stator base on differences of singularlyvalue[J]. Journal of Aerospace Power, 2009, 24(10)：2296-2301.

[11]李月, 杨宝俊. 混沌振子系统(L-Y)与检测[M]. 北京:科学出版社, 2007:12.

[12]Shin K, Hammond J K, White P R. Iterative SVD method for noise reduction of low-dimensional chaotic time series[J]. Mechanical Systems and Signal Processing, 1999,13(1): 115-124.

[13]徐锋, 刘云飞. 基于中值滤波-奇异值分解的胶合板拉伸声发射信号降噪方法研究[J]. 振动与冲击, 2011, 30(12): 135-140.

XU Feng, LIU Yun-fei. Noise reduction of acoustic emission signal in aplaywood based on median filtering-singular value decomposition[J]. Journal of Vibration and Shock, 2011, 30(12):135-140.

[14]赵学智, 向可, 叶邦彦,等. 基于二次样条小波细节信号峰值的有效奇异值确定[J]. 振动与冲击，2010, 29(11): 6-12.

ZHAO Xue-zhi, XIANG Ke, YE Bang-yan , et al. Determination of effective singular values based on detailed signal peak with quadratic spline wavelet [J]. Journal of Vibration and Shock, 2010, 29(11): 6-12.

[15]赵学智, 叶邦彦. 奇异值差分谱理论及其在车床主轴箱故障诊断中的应用[J]. 机械工程学报, 2010, 46(1):100-108.

ZHAO Xue-zhi, YE Bang-yan. Processing of milling force signal and isolation of state information of milling machine based on singular valuedecomposition[J]. Chinese Journal of Mechanical Engineering,2007, 43(06): 169-174.

[16]Shorter P J, Langley R S. Vibro-acoustic analysis of complex systems[J].Journal of Sound and Vibration, 2005, 288:669-699.