基于时频分析的运载火箭故障特征提取技术

崔桂磊 唐永华

(宇航动力学国家重点实验室，西安 710043)

基于时频分析的运载火箭故障特征提取技术

崔桂磊 唐永华

(宇航动力学国家重点实验室，西安 710043)

为尽可能从运载火箭飞行试验数据中获取更多有价值的信息，提出了新的火箭振动参数故障特征提取模型，用于火箭振动问题分析.构建了三维傅里叶谱，用于高效地分析全程观测信号频谱的时变趋势，研究了Wigner-Ville分布的标准算法与快速算法，用于精确刻画振动参数感兴趣数据段的频谱随时间的分布情况.试验结果表明:模型及方法处理火箭遥测数据直观、高效、精确，可以应用于实测数据分析.研究时频分析方法对试验数据进行更深层次的分析处理，对设计高可靠、高性能的运载火箭，具有重要意义.

运载火箭;遥测数据;故障定位;时频分析;特征提取

飞行试验遥测数据反映运载火箭的工作环境和各系统、各部件的工作状态，对于分析飞行试验过程是否正常以及分析飞行中出现的各种异常现象有着重要作用［1］.目前，火箭故障诊断的研究主要包括四个方面的内容:故障机理研究、状态信号采集、信号分析处理和故障特征提取、故障状态识别和诊断决策，其中的核心问题是故障信号的特征提取，其处理结果是故障诊断的重要依据.无论从地面试车状态还是飞行的实际来看，火箭遥测数据中往往含有大量的时变、短时突发性质的成分，其频率在故障发生时刻是跳变的，这使得目前基于传统线性和平稳理论发展起来的遥测数据处理技术如傅里叶变换、现代谱估计等［2］，难以有效地揭示故障信号的非平稳特征，因为它们都是对测量数据的整体变换，即对信号的表征或是完全在时域，或是完全在频域，只能对信号进行时域或者频域分析，作为频域表示的功率谱并不能说明其中某种频率分量出现在什么时候及其变化情况，故无法满足遥测故障信号的分析处理需要.在这种情况下，必须研究建立新的故障特征提取模型进行飞行试验数据的处理，以提高故障信号特征提取的品质，为进一步提高对试验数据的认识水平，进行型号任务火箭状态分析与故障定位提供保障.

1 两种时频分析方法

时频分析被认为是近年来对以傅里叶变换为基础的线性和平稳数据分析的一个重大突破，产生的原因主要是由于傅里叶变换在分析非平稳信号方面存在局限性，不能有效分析信号的瞬时特性［3－4］.时频分析提供了信号从时域到时间-频率域的变换，通过二维和三维时频分布图，能够在时频平面上直观表示出故障信号中各个分量的时间关联谱特性，在每一时刻指示出信号在瞬时频率附近的能量聚集情况，从而同时掌握信号的时域、频域信息，清楚地了解信号的频谱、能量谱是如何随时间变化的［5－6］.作为一种高效、高分辨率的数据处理技术，时频分析理应在飞行试验数据的使用和挖掘方面得到更好的运用.时频分析技术的基本思想是设计时间和频率联合函数，将一维时域信号映射到二维时频平面，在二维时间-频域平面全面分析观测信号的时频联合特征，实现对信号频谱或能量谱的时间定位［7－9］.按照时频联合函数的不同可以分为线性时频分析和非线性时频分析两类，其时频分辨率和聚集性各不相同［10－11］.因此工程实际中要研究最优方法来处理遥测信号，以获得高效、高分辨率的火箭飞行环境分析效果.

1.1 三维傅里叶谱

结合运载火箭遥测振动参数的特性，提出了全程数据三维傅里叶谱处理方案:

1)将振动参数在飞行试验的全程时域观测数据以1 024点为一个数据帧，划分为连续的若干帧信号，每个数据帧在时域占用的时间片长度为Δt=1024/fsr，fsr为该参数的采样率.

2)对每个数据帧进行快速傅里叶变换FFT(Fast Fourier Transform)变换得到其频谱结构，频谱的分辨率Δf=fsr/1024.

3)构建三维傅里叶谱的原理如图1所示，图中X为时间序列，以起测时间T0为时间原点，以时间片长度Δt为时间间隔逐帧增加;Y为频率序列，以谱分析的下限频率(一般为0到5Hz)为原点，以频谱分辨率Δf为频率间隔;Z为频谱幅值.

1.2 Wigner-Ville分布

在所有非线性时频表示中，Wigner-Ville分布(WVD)［12］是非常重要的一类二次时频分布，WVD的时频带宽积达到了Heisenberg不确定性原理给出的下界，使得WVD比线性时频分析具有更高的时频分辨率，时频聚集性比较好，并在许多方面如地震信号处理、故障特征提取、高光谱探测等［13－15］取得了成功的应用.

图1 三维傅里叶谱构造原理

1.2.1 WVD 的原理

设观测信号为s(t)，它的WVD表示为

式中，t为时间变量;f为频率变量;x(t)是信号s(t)的解析信号，即 x(t)=s(t)+j H［s(t)］，其中H［］表示希尔伯特变换;*表示复共轭.

式(1)中，信号x(t)出现了两次，即WVD以双线性或非线性的形式将信号从时域转换到时频域，描述了信号能量谱随时间、频率的变化规律，且不含任何窗函数，这样就避免了线性时频分析中时间和频率分辨率的相互牵制.

1.2.2 WVD标准算法与快速算法

1)WVD标准算法.

设观测信号为 s(n)，n=1，2，…，N，其离散解析信号为 x(n)，n=1，2，…，N，则其离散 WVD 可以表示为

式中，WM=e－j2π/M，M 为窗口长度.

分析式(2)可知，离散WVD可以通过FFT变换来实现.

式(2)可以视为WVD标准算法.标准算法将长度为N的观测信号作为一个整体进行WVD变换.作为信号的二次变换，WVD具有分辨率高、时频聚集性好的优点，但随着N的增大，时效性变差.

2)WVD快速算法.

当遥测信号长度N较大时，应用WVD标准算法进行时频分析，其处理效率往往不能满足工程快速分析的要求，因此实际中提出了快速算法方案.具体方案如下:①将待处理数据x(n)按1024点为1帧划分为n=N/1024帧，xi(n)为第i帧数据;②对划分好的数据逐帧进行WVD变换Φ［xi(n)］=Φxi(n，k);③在时频域将每帧处理结果进行拼接.

快速算法可以有效提高长数据序列的处理效率，与之对应的是时频分辨率有所降低，且频谱在帧间分割处会出现断痕，但一般不会影响对振动参数的正确分析.

2 火箭振动参数故障特征提取模型

根据测量目的不同，火箭上多个部位安置了不同类型的力学环境参数传感器，一般可分为振动参数、冲击参数、噪声参数和脉动压力参数等，其中振动参数又包括POGO振动、低频振动和高频振动［1］.由于运载火箭发射成本相当高昂，不可能过于频繁地进行试验，因此遥测数据具有很高的价值.环境参数处理的任务就是对测量数据进行分析，描述它们在时域、频域的基本特性，为火箭的设计评测及故障诊断提供可靠的依据.

对于长度为N的火箭遥测数据，经过大量数据试验，建立了如图2所示的振动参数故障特征提取模型，描述如下:

1)全程数据预处理.读取待分析参数的全程数据文件，对全程数据进行中心化、滤波等预处理工作，并绘制其时域波形图.

图2 火箭振动参数故障特征提取模型

2)全程数据三维傅里叶谱分析.全程数据以每1024点为1个数据帧，逐帧进行FFT变换，并在频率-时间-幅值三维空间绘制全程三维傅里叶谱，观测全程数据频谱变化的大概趋势.

3)数据选段.根据遥测数据处理要求及全程三维傅里叶谱，选取感兴趣的数据段落，长度为1024×m，m为帧数.对选段数据进行傅里叶变换，观察分析频谱结构，根据频谱信息进行时域数据滤波.

4)选段数据的WVD分析处理.根据帧数m是否小于等于4，选取WVD标准算法或快速算法.观察WVD二维时频谱、三维能量谱，分析故障信号特征频率随时间的分布情况.

5)进一步分析处理.若需要获取更细微的信息，则回到步骤3)，进一步进行数据选段并进行WVD时频分析.

3 火箭振动参数故障特征提取实例

以某次试验任务某型号运载火箭低频振动参数ZD01观测数据为例，分析如下:

1)全程数据预处理.遥测参数ZD01测量时间段为 －29.10 ～200.03 s，全程数据经中心化及5～100Hz滤波后，其时域波形如图3所示.

图3 振动参数全程数据时域图

图4 振动参数全程三维傅里叶谱

2)全程数据三维傅里叶谱分析.通过图4所示的全程三维傅里叶谱，可观察到120 s及150 s附近频率成分极为丰富.同时可以大致判断参数ZD01在0～100 s有较为明显的40Hz左右的振动频率存在，这就是型号任务专家高度重视的低频振动问题.但由于该方法为逐帧分析绘图，每帧1024点，因此时间分辨率与频率分辨率均较低，仅适用于遥测全程数据谱结构的粗略分析.

3)数据选段.选取参数 ZD01在 －1.23～101.15 s范围内的N=32768点(32帧)遥测数据进行时频分析处理，其时域波形如图5所示.传统傅里叶变换的频谱如图6所示，可以看出存在明显的40Hz频率信号，但是，图6的傅里叶频谱是对－1.23～101.15 s数据从时域到频域的整体变换，它只能说明在 －1.23 ～101.15 s中40Hz信号存在，却不能进一步给出关于40 Hz频率出现的时间以及40Hz频谱、能量谱随时间的分布情况.

图5 选段数据时域图

图6 选段数据频谱图

4)选段数据的WVD分析处理.对选段数据做WVD变换，将信号分析从时域转换到时间-频率域，在时频域全面观测信号的联合时频特性.由于选取的数据较长，因此应用WVD快速算法进行时频分析处理.图7为选段数据WVD的二维时频谱图，可以清晰地看出，参数ZD01在10～60 s存在40Hz信号;图8为选段数据WVD的三维能量谱图，进一步描述了40 Hz频率信号的能量随时间的分布情况.

5)进一步分析处理.选取参数ZD01在10～60 s范围内的遥测数据可以获得40 Hz频率更精确的时变特性;选取120 s与150 s附近的遥测数据，可以分别得到火箭一级关机与一二级分离时刻的振动特性，本文不再赘述.

6)振动参数ZD01的时变特性分析结论.由全程三维傅里叶谱、WVD的二维时频谱及三维能量谱全面观测信号的联合时频特性，可以准确判断出40Hz振动频率发生与消失的时间及其变化趋势.

图7 选段数据WVD的二维时频谱

图8 选段数据WVD的三维能量谱

4 结论

1)着眼于长征系列运载火箭振动问题分析的需要，建立了新的故障特征提取模型，从联合时频分布角度直观描述了振动参数特征频率的频谱与能量谱随时间的变换趋势，提高了对试验数据的认识水平.

2)针对遥测振动参数时变特性，构建了三维傅里叶谱，研究了Wigner-Ville分布标准算法及快速算法，实现了箭遥振动参数从单纯的时域或频域分析到时频域分析的转变，提高了数据处理能力.

3)实测数据试验表明，新处理模型直观、高效，可提供更精确、更丰富的诊断信息，其方法结果参考价值更高，能够为型号部门进行故障状态识别和诊断决策提供更为科学的依据.

References)

［1］陈以恩，李连周，唐永华，等.遥测数据处理［M］.北京:国防工业出版社，2002:3－20 Chen Yien，Li Lianzhou，Tang Yonghua，et al.Telemetry data processing［M］.Beijing:National Defence Industry Press，2002:3－20(in Chinese)

［2］GJB2238A—2004遥测数据处理［S］GJB2238A—2004 Telemetry data processing［S］(in Chinese)

［3］ Chaparro L F，Akan A.Special issue on time-frequency signal analysis and its applications［J］.Journal of the Franklin Institute，2000，337(4):287 －289

［4］ Samantaray L，Dash M，Panda R.A review on time-frequency，time-scale and scale-frequency domain signal analysis［J］.IETE Journal of Research，2005，51(4):287 －293

［5］ Stankovic S.Time-frequency analysis and its application in digital watermarking［J］.EURASIP Journal on Advances in Signal Processing，2010(2):1 －20

［6］ Sejdic E，Djurovic I，Jiang J.Time-frequency feature representation using energy concentration:an overview of recent advances［J］.Digital Signal Processing，2009，19(1):153 － 183

［7］ Ristic B，Boashash B.Kernel design for time-frequency signal analysis using the Radon transform［J］.IEEE Transactionson Signal Processing，1993，41(5):1996 －2008

［8］ Barkat B，Boashash B.A high-resolution quadratic time-frequency distribution for multicomponent signals analysis［J］.IEEE Transactions on Signal Processing，2001，49(10):2232 －2239

［9］ Shafi，I，Ahmad J，Shah，S I，et al.Techniques to obtain good resolution and concentrated time-frequency distributions:a review［J］.EURASIP Journal on Advances in Signal Processing，2009(2):1－43

［10］ Hlawatsch F，Boudreaux-Bartels G F.Linear and quadratic timefrequency signal representations［J］.IEEE Signal Processing，1992，9(2):21 －67

［11］ Jeffrey C，Patrick F.Virtues and vices of quadratic time-frequency distributions［J］.IEEE Transactions on Signal Processing，2000，48(9):2641 －2650

［12］ Cohen L.Time-frequency distributions-a review［J］.Proceedings of the IEEE，1989，77(7):941 －981

［13］ Boashash B，Black P J.An efficient real-time implementation of the Wigner-Ville distribution［J］.IEEE Transactions on Acoustics Speech and Signal Processing，1987，35(11):1611 －1618

［14］ Xu Y N，Zhao Y，Jin C F，et al.Salient target detection based on pseudo-Wigner-Ville distribution and Renyi entropy［J］.Optics Letters，2010，35(4):475 －477

［15］ Climente-Alarcon V，Antonino-Daviu JA，Riera M，et al.Wigner-Ville distribution for the detection of high-order harmonics due to rotor asymmetries［C］//IEEE International Symposium on Diagnostics for Electric Machines，Power Electronics and Drives.［s.l.］:IEEE，2009:119 －124

(编 辑:李 晶)

Rocket fault feature extraction based on time-frequency analysis

Cui Guilei Tang Yonghua

(State Key Laboratory of Astronautic Dynamics，Xian 710043，China)

To obtain more valuable information from rocket flight test data asmuch as possible，a fault feature extraction model was proposed to the need of rocket vibration analysis.3D Fourier spectrum was constucted to analyze the time-varying trend of entire observation signal efficiently.Standard and fast algorithm of Wigner-Ville distribution was studied to describe the frequency spectrum of interesting data segments accurately.Themodel and methodswere illustrated on rocket telemetry data.The results prove the model's direct-vision，high efficiency，and accuracy.It's of great significance for designing a rocket with high reliability and performance to analyze the flight test data more deeply，based on time-frequency analysis methods.

rocket;telemetry data;fault location;time-frequency analysis;feature extraction

V 557+.3;TJ 760.6+28;TN 911.72

A

1001-5965(2012)06-0741-04

2011-03-10;网络出版时间:2012-06-15 15:42

www.cnki.net/kcms/detail/11.2625.V.20120615.1542.016.htm l

总装预研项目(2009SY2603009)

崔桂磊(1976 －)，男，山东青岛人，工程师，CGL11@mails.tsinghua.edu.cn.