正弦掃頻振動試驗結(jié)構(gòu)響應(yīng)的高精度頻譜分析方法

王營營，白 楊，徐曉輝，王繼虎，李霖圣

（1.上海衛(wèi)星裝備研究所，上海200240；2.上海航天技術(shù)研究院，上海201109）

0 引言

正弦掃頻振動試驗是航天器研制過程中常開展的力學(xué)環(huán)境模擬試驗，其試驗結(jié)果對于評估航天器耐受振動環(huán)境的能力、改進航天器在軌可靠性設(shè)計等起著重要作用。

力學(xué)環(huán)境試驗數(shù)據(jù)采集和分析處理系統(tǒng)具有對應(yīng)試驗時間短、測量通道多、數(shù)據(jù)流量大、實時性要求高、可靠性要求高等特點。由于單套數(shù)據(jù)采集系統(tǒng)的通道數(shù)發(fā)展受到硬件性能的限制，以分布集成式系統(tǒng)代替超大系統(tǒng)模式成為多通道大型數(shù)據(jù)采集系統(tǒng)的發(fā)展趨勢[1]，而大量分布式采集和存儲的結(jié)構(gòu)響應(yīng)時域數(shù)據(jù)對數(shù)據(jù)后處理能力提出了較高要求。

航天器上安裝的力學(xué)環(huán)境測量系統(tǒng)可以測量真實的主動段動力學(xué)環(huán)境和在軌微振動力學(xué)環(huán)境，對航天器的力學(xué)環(huán)境優(yōu)化和設(shè)計有重要參考價值[2]。地面振動試驗期間進行測試，并和地面數(shù)據(jù)采集系統(tǒng)得到的幅值譜進行對比分析，可以驗證在軌力學(xué)環(huán)境測量系統(tǒng)的有效性[3]。由于在軌力學(xué)環(huán)境測量系統(tǒng)在航天器上獨立工作，和地面振動控制儀不能實現(xiàn)同步，對其采集的時域振動響應(yīng)數(shù)據(jù)進行準(zhǔn)確的頻譜分析才能實現(xiàn)數(shù)據(jù)對比和開展系統(tǒng)有效性研究。

地面振動試驗數(shù)據(jù)采集系統(tǒng)通常利用振動控制儀輸出的COLA（constant output level adapter）信號進行實時頻率辨識，從而實現(xiàn)同步數(shù)據(jù)采集和在線數(shù)據(jù)處理。實時處理振動響應(yīng)數(shù)據(jù)對系統(tǒng)硬件要求高，而系統(tǒng)的通道數(shù)發(fā)展受到硬件制約，因此不依賴COLA 信號的正弦振動試驗時域數(shù)據(jù)分析方法在結(jié)構(gòu)響應(yīng)分析中有重大研究意義。

有效的正弦振動響應(yīng)信號分析基于對適用于非平穩(wěn)信號的時頻分析方法的探索。短時傅里葉變換（STFT）定義的譜圖是簡單、直觀的時頻表示，但其受到測不準(zhǔn)原理的制約，不能同時兼顧時間分辨率和頻率分辨率[4]。使用離散傅里葉變換（discrete Fourier transform,DFT）分析振動響應(yīng)信號時，時域信號的整周期截斷和同步采樣通常難以實現(xiàn)，頻譜泄漏會造成嚴重的分析誤差[5]。一般利用窗函數(shù)和插值方法來減少頻譜泄漏[6]，通過調(diào)整窗寬和采樣率等技術(shù)可以在一定程度上減小誤差[7]，但是窗函數(shù)的引入會對幅值分析的精度產(chǎn)生影響[8]。針對DFT方法中非同步采樣信號的幅值和頻率分析精度不高的問題[9]，對泄漏頻譜進行插值校正[10]或推導(dǎo)出真實頻譜[11]都可以達到減小誤差的效果。針對信號的非整周期截斷問題，采用殘周期正弦擬合方法[12]可以獲取低頻信號的振動參數(shù)?；贖ilbert 變換和濾波算法可以利用驅(qū)動信號[13]或構(gòu)造參考信號[14]分析得到頻譜曲線，但是濾波處理的精度和結(jié)構(gòu)特性有關(guān)，并且分段濾波會導(dǎo)致局部的幅值誤差較大。單幀F(xiàn)FT 法[15]可以有效識別結(jié)構(gòu)動特性和頻響曲線，但是無法得到準(zhǔn)確的幅值譜?；谧赃m應(yīng)帶通濾波的STFT方法[16]能有效改善STFT幅值分析的精度，但與地面數(shù)據(jù)采集系統(tǒng)的處理結(jié)果相比有較大誤差。

本文針對不依賴COLA 信號的正弦振動響應(yīng)分析，提出一種基于掃頻規(guī)律的頻率跟蹤DFT 方法，根據(jù)正弦掃頻試驗的頻率變化特性[17]，基于過零點檢測和最小二乘法進行信號頻率辨識，同時在較小的時間尺度上將正弦掃頻信號近似為平穩(wěn)周期信號，根據(jù)基頻對時域振動響應(yīng)信號進行整周期截斷和DFT頻譜分析，旨在提高正弦振動響應(yīng)頻率和幅值的分析精度。

1 時頻分析理論

1.1 正弦掃頻振動試驗概述

正弦掃頻振動試驗是按照規(guī)定的振動量級，在一定頻率范圍內(nèi)連續(xù)改變振動頻率來激勵試驗件。其振動頻率的變化率稱為掃描速率。

試驗激振頻率f、試驗時間t和掃描速率β的關(guān)系在線性掃頻中表示為

在對數(shù)掃頻中表示為

將振動試驗的起始頻率記為f0，起始時間記為t0，利用式(2)可得到對數(shù)掃頻試驗的實時頻率f和時間t的關(guān)系為

若取β=k×ln 2，則可得到工程中常用的倍頻程表示形式

式中k為倍頻程掃描速率。

激振頻率f與時間t在線性掃頻和對數(shù)掃頻試驗中均有比較簡單的函數(shù)關(guān)系，對于復(fù)雜的頻率掃描方法，可以用曲線函數(shù)、分段函數(shù)或者離散點序列等方式來描述。正弦掃頻試驗中的頻率掃描規(guī)律是已知的，只有起始時間t0是由數(shù)據(jù)采集器定義的，只需確定t0就能按照頻率掃描規(guī)律計算出任意時刻的激振頻率。

1.2 頻譜分析方法

正弦振動試驗地面數(shù)據(jù)采集系統(tǒng)通過COLA信號得到振動臺的實時激振頻率，系統(tǒng)硬件以高采樣率進行時域數(shù)據(jù)采集、實時分析處理、數(shù)據(jù)傳輸和存儲。振動響應(yīng)幅值譜的計算方法主要有絕對峰值法、平均值峰值法、有效值峰值法和濾波方式峰值法。除絕對峰值法外，其他3種方法都有一定的濾波效果；濾波方式峰值法的分析速度快，處理結(jié)果曲線平滑，并且可以獲得相位信息，是工程中常用的分析方法。

濾波方式峰值法將正弦振動信號表示為

式中：N為數(shù)據(jù)塊中的數(shù)據(jù)個數(shù)；Δt為采樣周期。

地面數(shù)據(jù)采集系統(tǒng)的硬件采樣率通常遠高于信號基頻：以比利時LMS公司的SCADASⅢ系統(tǒng)為例，其前端單通道采樣率可達到204.8 kHz，而振動試驗頻率上限一般小于2 kHz。系統(tǒng)在1個信號周期內(nèi)即可采集到足夠的數(shù)據(jù)，然后通過最小二乘法和式(7)得到ac和as的精確估計值，從而計算出信號幅值A(chǔ)(ω1)。

使用振動響應(yīng)時域數(shù)據(jù)后處理的方法進行頻譜分析主要有2個難點：1）通常沒有COLA 信號作為參考，難以進行精確的頻率辨識；2）受時域數(shù)據(jù)傳輸和存儲過程的硬件性能限制，實際的采樣率較低，通過式(7)進行幅值估計時，分析結(jié)果的抗干擾能力較差。

2 頻率跟蹤DFT方法

2.1 正弦掃頻信號的頻率識別

通常正弦振動試驗激勵頻率和結(jié)構(gòu)響應(yīng)頻率的相關(guān)性很高，因此，由過零點檢測方法得到的數(shù)據(jù)點(tr,lnfr)應(yīng)該落在式(8)所描述的直線附近。如圖1所示，對過零點檢測法計算出的數(shù)據(jù)點序列使用最小二乘法即可得到lnf和時間t的直線表達式；利用式(8)和起始頻率f0可得到起始時間t0，然后根據(jù)式(4)進行頻譜分析中的頻率辨識。

對于試驗頻率單調(diào)遞增或單調(diào)遞減的掃頻方式，只需要1對精確的時間?頻率值就能計算出起始時間t0，但考慮由過零點檢測方法得到的頻率值容易受到干擾信號影響，故需要取響應(yīng)信號的多個過零點檢測值并利用最小二乘法減小分析誤差。試驗頻率不符合單調(diào)遞增或單調(diào)遞減規(guī)律時，可以先估計起始時間t0的大致范圍，使得t0的取值范圍小于試驗頻率的變化周期。

結(jié)構(gòu)的固有頻率特性和信號干擾對過零點檢測方法的頻率計算精度影響較大，一般可以通過濾波和設(shè)置閾值的方法排除數(shù)據(jù)畸變點。受到非同步采樣的影響，在采樣率較低時，高頻信號的頻率計算誤差也比較大，在最小二乘法中只選取低頻段的數(shù)據(jù)點可有效減小誤差，提高頻率辨識精度。

2.2 頻率識別算例

對于二次函數(shù)掃頻正弦振動，瞬時頻率fi(t)是時間t的二次函數(shù)。例如，式

的含義為：0 s—30 s的時間內(nèi)，瞬時頻率從100 Hz按照二次函數(shù)規(guī)律連續(xù)變化到1 Hz；在30 s—60 s的時間內(nèi)，再從1 Hz 連續(xù)變化到100 Hz。

利用式(12)的掃頻規(guī)律可以構(gòu)造加速度幅值為1g的掃頻數(shù)據(jù)。為模擬數(shù)據(jù)采集器時間和掃頻時間不能同步的情況，在時域掃頻數(shù)據(jù)前后分別增加5 s 的定頻振動（加速度幅值1g、頻率100 Hz）。將掃頻起始時間記為t0，顯然在掃頻數(shù)據(jù)前面增加了5 s的定頻數(shù)據(jù)后，有t0=5。在式(12)中引入起始時間t0，將其改寫為

下面假設(shè)t0是未知量，并使用過零點檢測法和最小二乘法對構(gòu)造的時域掃頻數(shù)據(jù)進行分析，來詳細說明求解t0的過程。在開始求解t0之前，將圖2所示的沖擊響應(yīng)信號疊加到時域掃頻數(shù)據(jù)中，為過零點檢測方法引入干擾項，最終得到的時域掃頻數(shù)據(jù)見圖3，時間長度70 s，STFT得到的頻譜如圖4 所示。

第1步，在時域數(shù)據(jù)中每隔2 s（分別從時刻0 s、2 s、4 s、 ···、68 s開始）使用過零點檢測方法按照式(9)向后尋找2個過零點，2個過零點時間差的倒數(shù)為頻率值fi，2個過零點時刻的平均值為時刻ti。在70 s的掃頻數(shù)據(jù)總時長中，可以計算得到35組時間?頻率值(ti,fi)，其中1≤i≤35。

圖2 沖擊響應(yīng)信號Fig.2 Shock response signal

圖3 二次函數(shù)掃頻振動時域響應(yīng)Fig.3 The time domain response of quadratic sine-swept vibration

圖4 二次函數(shù)掃頻振動頻譜Fig.4 The spectrum of quadratic sine-swept vibration

第2步，選取低頻數(shù)據(jù)計算均方根誤差。由于低頻部分每個信號周期內(nèi)的采樣點數(shù)目較多，過零點分析所得頻率值fi的誤差小，所以優(yōu)先使用低頻數(shù)據(jù)。本文取1≤f≤5，根據(jù)式(13)可得ti需要滿足條件：由于0≤t≤70，根據(jù)式(13)可得起始時間t0可能的取值范圍為[0,10]，在該范圍內(nèi)取t0為任意值，通過式(14)即可選取出[ta,tb]范圍內(nèi)的過零點分析值(ti,fi)，其中m≤i≤n。圖5所示為：t0=1的情況下m=14、n=18；圖6所示為：t0=9的情況下m=18、n=23。

圖5 瞬時頻率曲線（t0=1 s）Fig.5 Thecurve of instantaneousfrequency (t0=1 s)

圖6 瞬時頻率曲線（t0=9 s）Fig.6 Thecurve of instantaneousfrequency (t0=9 s)

的過零點分析值(ti,fi)剔除；將區(qū)間[ta,tb]內(nèi)剩余的N個過零點分析值記為(tj,fj)，其中1≤j≤N，重新計算均方根誤差

對于任意的t0值，均可按照上述方法計算出相應(yīng)的均方根誤差ERMS。

利用最小二乘法的原理，求解t0最優(yōu)解的問題可以描述為：在[0,10]范圍內(nèi)求t0的值，使得均方根誤差ERMS的值最小。求解該問題最直接的方法是在[0,10]范圍內(nèi)以采樣周期為步長求解所有ERMS，通過尋找ERMS的最小值來確定對應(yīng)的t0。該問題也可以使用工程數(shù)值優(yōu)化算法快速求得最優(yōu)解（如圖7所示），其值為t0=5。

圖7 最小二乘法計算結(jié)果Fig.7 The resultsof least squaresmethod

2.3 時頻分析方法

利用式(21)可得到ac和as的值，再根據(jù)式(6)即可得到該頻率處的信號幅值。該方法和濾波方式峰值法相似，在計算幅值的同時可獲取相位信息。

上述方法基于響應(yīng)信號的頻率辨識和頻率跟蹤，并使用整周期截斷有效減小了DFT頻譜分析誤差。由于數(shù)據(jù)截斷時將信號頻率視為常數(shù)，所以頻率變化太快或數(shù)據(jù)長度N較大時會產(chǎn)生截斷誤差。另一方面，由于信號的非同步采樣導(dǎo)致譜線頻率fk不會精確等于信號基頻，所以采樣頻率相對信號基頻太小或數(shù)據(jù)長度N較小時會增大頻譜分析誤差。航天器正弦振動試驗中，為了激發(fā)結(jié)構(gòu)在試驗頻率的穩(wěn)態(tài)響應(yīng)，掃描速率通常不宜太快，取周期數(shù)M=4可以有效地控制誤差。

3 對比分析及驗證

3.1 與構(gòu)造的正弦掃頻試驗數(shù)據(jù)幅值譜理論值的對比

算法構(gòu)造一個正弦掃頻數(shù)據(jù)序列，使用頻率跟蹤DFT 方法進行時頻分析，并將計算結(jié)果和理論值進行對比，驗證第2章所述方法的計算精度。

構(gòu)造正弦掃頻試驗數(shù)據(jù)，其采樣率為1000 Hz，頻率范圍1～200 Hz，掃描速率4 oct/min，加速度幅值1g。用頻率跟蹤DFT 方法進行時頻分析：頻率范圍2～100 Hz，頻率分辨率為每倍頻程100譜線。對構(gòu)造的時域數(shù)據(jù)進行過零點檢測并選取2～5 Hz 的頻率范圍進行頻率識別，頻域幅值譜的理論值和分析結(jié)果的曲線見圖8，具體結(jié)果對比和誤差數(shù)據(jù)見表1。

由圖8可見，幅值譜的高頻部分誤差稍大，這主要是因為采樣頻率相對高頻段信號的基頻較小，導(dǎo)致頻譜分析中非同步采樣誤差的影響較大。由表1可見，頻率跟蹤DFT 方法的分析結(jié)果和理論值間的誤差很小，幅值相對誤差最大為0.69%。

圖8 幅值譜理論值和頻率跟蹤DFT 方法結(jié)果Fig.8 The theoretical value of amplitude spectrum and that obtained by frequency tracking DFT

表1 頻率跟蹤DFT方法處理結(jié)果和理論值對比Table 1 Comparison between results obtained by frequency tracking DFT and theoretical value

采樣率為400 Hz 時，對掃描速率為4 oct/min 的對數(shù)掃頻和式(13)描述的二次函數(shù)掃頻分別使用本文方法進行分析，幅值分析結(jié)果與理論值之間的誤差見圖9（其中，左圖為頻域幅值譜的理論值和分析結(jié)果的曲線，右圖為相對誤差數(shù)據(jù)）?？梢钥吹剑赫`差大小和1個信號周期內(nèi)的采樣點個數(shù)（即采樣頻率與信號頻率的比值）相關(guān)——采樣頻率大于掃頻頻率的4倍時分析誤差在5%以內(nèi)，采樣頻率大于掃頻頻率的10倍時分析誤差在1.5%以內(nèi)。

圖9 頻率跟蹤DFT 方法所得幅值譜的誤差Fig.9 The errorsof theamplitude spectrum obtained by frequency tracking DFT ascompared with the theoretical values

3.2 與地面數(shù)據(jù)采集系統(tǒng)幅值譜處理結(jié)果對比

使用頻率跟蹤DFT 方法分析某航天器的正弦振動試驗數(shù)據(jù)，并和地面數(shù)據(jù)采集系統(tǒng)的濾波方式峰值法處理結(jié)果進行對比。

結(jié)構(gòu)振動響應(yīng)測量點有3個：1號測點在航天器上部，2號測點在航天器中部，3號測點在航天器底部。試驗的頻率范圍為5～100 Hz，掃描速率為4 oct/min。振動激勵方向為航天器x方向，同時采集了x向和y向的振動響應(yīng)時域數(shù)據(jù)和頻域幅值譜數(shù)據(jù)。時域數(shù)據(jù)采樣率為400 Hz，加速度響應(yīng)曲線如圖10所示（其中，1x、1y分別代表1號測點的x向、y向振動響應(yīng)，其他依此類推）。頻率分辨率為每倍頻程100 譜線。

圖10 某航天器振動響應(yīng)時域數(shù)據(jù)Fig.10 The time domain vibration response of the spacecraft

由于x向是振動激勵方向，所以對x向響應(yīng)數(shù)據(jù)進行頻率識別的精度較高。對3x測點響應(yīng)數(shù)據(jù)的頻率辨識結(jié)果如圖11所示，可得到掃描起始時間為16.865 s。振動響應(yīng)數(shù)據(jù)采集時使用的是同一套數(shù)據(jù)采集系統(tǒng)，因此各通道的時域數(shù)據(jù)對應(yīng)的時間是同步的，由3x測點得到的掃描起始時間適用于所有通道。

圖11 頻率辨識結(jié)果Fig.11 The frequency identification results

由圖11可見，過零點檢測法得到的頻率值在35～40 s之間存在明顯的異常數(shù)據(jù)，這和結(jié)構(gòu)的一階頻率特性有關(guān)。另外，結(jié)構(gòu)響應(yīng)中的撞擊信號也會使該頻率檢測方法計算得到異常值。使用低頻段（5～7 Hz）范圍的數(shù)據(jù)進行最小二乘法頻率辨識，可有效避免響應(yīng)信號中的異常干擾，顯著提高頻率識別精度。

根據(jù)振動試驗的掃描起始時間進行頻率識別，使用頻率跟蹤DFT 方法對振動響應(yīng)數(shù)據(jù)進行頻譜分析，所得頻譜曲線見圖12，可見分析所得幅值譜和地面數(shù)據(jù)采集系統(tǒng)的處理結(jié)果基本一致，具體誤差數(shù)據(jù)見表2。表2列出了幅值譜曲線的最大峰值及其對應(yīng)的頻率值，對比顯示，頻率跟蹤DFT方法和地面數(shù)據(jù)采集系統(tǒng)所得結(jié)果之間的峰值誤差最大為2y測點的1.149 2%，振動激勵方向的峰值誤差最大為1x測點的-0.369 3%；峰值頻率誤差最大為-0.688 5%；幅值譜平均誤差最大為-0.000 7g，均方根誤差最大為0.035 3g，相關(guān)系數(shù)均在0.99以上。頻率跟蹤DFT 方法顯著減小了非整周期截斷引入的頻譜分析誤差，幅值分析精度和地面數(shù)據(jù)采集系統(tǒng)濾波方式峰值法的處理精度相差不大。同時，該方法不需要使用窗函數(shù)并且可以獲得準(zhǔn)確的幅值譜，與傳統(tǒng)的時頻分析方法相比具有明顯優(yōu)勢。

圖12 頻率跟蹤DFT 方法和地面數(shù)據(jù)采集系統(tǒng)所得幅值譜Fig.12 The amplitude spectrum obtained by frequency tracking DFT and the ground data acquisition system

表2 頻率跟蹤DFT 方法和地面數(shù)據(jù)采集系統(tǒng)所得結(jié)果對比Table 2 Comparison of resultsobtained by frequency tracking DFTand the ground data acquisition system

當(dāng)時域數(shù)據(jù)的采樣率較高時，非同步采樣誤差的影響變小，頻率跟蹤DFT 方法計算得到的幅值譜精度將會更高。當(dāng)一套數(shù)據(jù)采集系統(tǒng)同時采集的結(jié)構(gòu)響應(yīng)數(shù)據(jù)中沒有激振方向的數(shù)據(jù)時，需要選取和振動控制相關(guān)性較好的響應(yīng)數(shù)據(jù)進行頻率辨識，并進行必要的濾波處理。

4 結(jié)束語

本文提出了頻率跟蹤DFT方法，應(yīng)用于正弦掃頻振動試驗數(shù)據(jù)處理可以得到高精度的幅值譜。該方法不依賴振動控制儀輸出的COLA 信號，也不需要數(shù)據(jù)采集系統(tǒng)對響應(yīng)信號進行在線頻譜分析，對系統(tǒng)的硬件性能要求比較低，有利于發(fā)展多通道分布式大型數(shù)據(jù)采集系統(tǒng)；在時域數(shù)據(jù)采樣率較低的情況下，能得到準(zhǔn)確的頻譜數(shù)據(jù)，分析結(jié)果和地面數(shù)據(jù)采集系統(tǒng)的處理結(jié)果間誤差很小，適用于振動環(huán)境試驗中航天器在軌力學(xué)環(huán)境測量系統(tǒng)的有效性驗證；使用整周期截斷技術(shù)避免了傳統(tǒng)時頻分析方法中的窗函數(shù)引入的分析誤差。如果對高頻段的振動響應(yīng)分析精度要求較高，可采用針對非同步采樣的頻譜校正算法進一步提高頻譜分析精度。值得一提的是，在結(jié)構(gòu)響應(yīng)信號失真嚴重的情況下，該方法的頻譜分析精度及適用性有待進一步研究與驗證。