黃繼堯,陳長興,凌云飛,藺向陽
(空軍工程大學基礎部,陜西西安 710043)
Duffing振子是一種用于信號處理的非線性混沌振子模型,它在混沌臨界狀態(tài)下對噪聲具有免疫性,但對初始條件極為敏感[1]??梢岳肈uffing振子的這些特性來檢測強噪聲背景下的微弱信號。自20世紀后期引入混沌振子系統(tǒng)以檢測微弱信號[2-4]以來,諸多學者在此領域進行了深入的研究。謝濤等證明了混沌振子檢測微弱信號的可靠性,為混沌檢測微弱信號提供了理論論據(jù)[5]。
李月等分析了不同類型噪聲對混沌系統(tǒng)的影響,證明了混沌系統(tǒng)對噪聲的免疫性[6];秦衛(wèi)陽等分析了Duffing方程的派生系統(tǒng)在同步過程中能夠?qū)崿F(xiàn)檢測2種信號中的微小差別[7];劉海波等提出一種利用單一振子實現(xiàn)以Duffing振子大周期態(tài)為基礎的測量頻率的方法[8]。常見的混沌狀態(tài)判定方法有Lyapunov指數(shù)法、龐加萊截面法、分維數(shù)計算法、Kolmogorov熵法、Melnikov法等[9]。與傳統(tǒng)的線性信號檢測方法相比,基于Duffing振子的微弱信號檢測方法可以實現(xiàn)非常低的信噪比(SNR),信噪比門限最低可達到-100 dB,并可得到傳統(tǒng)檢測方法無法得到的微弱待測信號,這為低信噪比的微弱信號檢測提供了一種新的方法[10]。
然而目前大多數(shù)基于Duffing振子的微弱信號檢測研究都是針對單頻諧波信號的檢測識別[11-14],而很少有研究關注微弱多頻信號的檢測。因此對基于Duffing振子的微弱多頻信號檢測進行深入研究是十分有必要的。本文提出了一種結(jié)合信號頻率截取預處理和Duffing振子檢測微弱多頻信號的方法。
Holmes型Duffing方程[15]一般形式為
(1)
式中:k為阻尼比;F為周期性策動力幅值;ω為角頻率;Fcosωt為周期性策動力。
當k保持不變時,隨著F的增大,存在一個臨界策動力振幅Fd。當F小于Fd時,系統(tǒng)處于混沌狀態(tài);當F大于Fd時,系統(tǒng)處于大周期態(tài)。在Duffing振子方程式(1)中,臨界值Fd≈0.825?;煦绾痛笾芷趹B(tài)的變化及策動力臨界值如圖1所示。
(a)混沌態(tài)
(b)周期態(tài)
(c)臨界點圖1 系統(tǒng)歷經(jīng)的狀態(tài)變化及分岔臨界圖
當Duffing系統(tǒng)處于從混沌狀態(tài)轉(zhuǎn)變?yōu)榇笾芷趹B(tài)的臨界狀態(tài)時,就統(tǒng)計特性而言,它不受高斯白噪聲的影響。由于對參數(shù) 的敏感性和對噪聲的免疫性,式(1)可以進一步表示為
(2)
式中:s(t)為待測信號;n(t)為噪聲;s(t)+n(t)為帶有噪聲的待測信號。
由式(2)可知,當式(1)中的F略小于臨界幅值Fd時,s(t)+n(t)輸入Duffing振子系統(tǒng),n(t)作為擾動項,則可以使用式(2)檢測待測信號。sn(t)=s(t)+n(t)是含噪待測信號。當s(t)包含角頻率為ω的信號分量時,系統(tǒng)將從混沌狀態(tài)轉(zhuǎn)換為大周期態(tài),這可以從系統(tǒng)響應的相平面軌跡獲得,然后便可以獲得待測信號的頻率。
傅里葉變換是一種分析信號的方法,可以把具有周期特征的非正弦信號分解為多個疊加的正弦信號,然后通過對時間連續(xù)信號的等間隔采樣,得到采樣頻率fs,并把采樣值轉(zhuǎn)換為數(shù)字序列。連續(xù)時間序列x(t)與其對應的連續(xù)非周期頻譜函數(shù)X(f)之間的轉(zhuǎn)換關系為:
(3)
當設定一個待測信號為hcos(ωt+φ),其采樣頻率為fs。ω為確定的角頻率,φ是初始相位。假設采樣信號向坐標軸左側(cè)移動dt=1/fs。根據(jù)公式T/dt=2π/dθ,就可以得到dθ=dt·2π/T=ω/fs,T為信號的周期,dθ為對應相位的增量。因此,將采樣信號向左挪動時間j等同于將待測信號由hcos(ωt+φ)變成hcos(ωt+φ+jdθ)。
在不考慮噪聲,即n(t)=0的情況下,設定待測信號尺度變換后(尺度變換系數(shù)R=ω)與策動力有相同的角頻率,則檢測待測信號的Duffing方程為
(4)
改變式(6)右側(cè)項的形式
(5)
從式(4)中,可以推出只有當(φ+jdθ)=2mπ(m是整數(shù))時,與待測信號和策動力相加的待測信號的振幅是最大的。
文獻[16]提出了一種基于Duffing振子的尺度變換微弱信號檢測模型,有效地解決了這些問題。檢測模型如式(6)所示:
(6)
這4個方程的輸入信號初始相位的檢測范圍能夠覆蓋整段頻率(-π,π),避免了單一方程只能檢測某一段初始相位的限制。
式(6)解決了基于Duffing振子微弱信號檢測的一些關鍵問題,然而目前僅對含噪聲的單頻諧波信號進行分析,而不涉及實際工程應用中常見的微弱多頻信號的檢測。因此本文將研究基于Duffing振子的微弱多頻信號檢測。
為了將式(6)中的檢測模型應用于微弱多頻信號,提出了一種信號頻率截取的預處理方法。其原理如下:根據(jù)Shannon采樣定理,試驗信號sn(t)的采樣頻率為fs,則其截止頻率為fs/2,為了便于后續(xù)分析,首先在頻域中截取待測信號。截取頻率的下限為fmin,上限為fmax,其滿足0 首先,利用快速傅里葉變換(FFT)將試驗信號sn(t)變換為頻域信號sn(f);其次,信號sn′(f)是由sn(f)通過頻率截取得到的,即離散譜數(shù)據(jù)的頻率范圍為(fmin,fmax);最后,將信號sn′(f)通過快速傅里葉逆變換(IFFT)轉(zhuǎn)變?yōu)闀r域信號sn′(t)。 上述頻率截獲過程可以表示為: (7) (8) (9) 頻率截取過程后,sn′(t)只包含滿足fmin 將信號頻率截獲方法與式(6)檢測模型相結(jié)合,針對微弱多頻信號,提出了一種新的待測信號提取模型。其過程具體如下: (1)對原始待測信號進行加噪處理,得到被噪聲污染的待測信號sn(t),采樣頻率為fs, 由強度為D的背景噪聲和微弱待測信號s1(t),s2(t),…,sn(t)構(gòu)成,微弱待測信號的角頻率和振幅分別為ω1,ω2,…,ωn和h1,h2,…,hn。 (2)確定原始待測信號的頻率測試范圍[ωmin,ωmax]。然后設定ω=ωmin,接著設置頻率增量dω,令ωi+1=ωi+dω。 (3)然后對待測信號sn(t)進行傅里葉變換后得到信號sn(f)再進行頻段為[ωmin,ωmax]的頻率截取預處理。新的信號在預處理后為sn′(f) ; (4)設定尺度系數(shù)Ri=ωi(即ω的值每增長一次,R的值都會相應的增長并與之相等),然后對信號進行頻率/時間尺度變換,之后得到在時間尺度t′下新的待測信號sn′(t)。變換后的采樣頻率為fsr=fs/Ri; (5)頻率截取預處理和尺度變換后,將新獲得的信號sn′(t)輸入到式(6)中。然后利用四階Runge-Kutta算法求解各方程,繪制各輸出響應的相平面軌跡圖。只要任何一個相平面軌跡處于大周期態(tài),則意味著待測信號sn(t)中的頻率分量ωk被檢測到,然后重復迭代過程,直到頻率滿足ω≥ωmax為止,過程如圖2所示。 圖2 信號檢測框圖 從而檢測到該頻率范圍內(nèi)的所有信號頻率分量,實現(xiàn)對未知多頻率信號的檢測。如果dω太小,檢測未知信號的計算量將大幅增加。選擇角頻率為ω=20 rad/s的待測信號進行測試,設定Δω為0.2 rad/s。選擇不同的Δω的倍數(shù)與計算時間的關系如圖3所示,選擇dω為2Δω有較快的運算速度,因此選擇dω為2Δω。 圖3 不同的Δω倍數(shù)與計算時間的關系圖 所以當檢測過程完成后,可以覆蓋完整的頻段[ωmin,ωmax]。迭代過程如圖4所示: 圖4 未知多頻信號檢測的迭代過程 設定待測信號的方程為 (10) 式中:ω1,ω2,ω3和h1,h2,h3分別為sn(t)中包含的3個頻率分量的角頻率和振幅;g(t)為均值為0,方差為1的高斯白噪聲。 設定待測信號的初始相位為0,待測信號的采樣頻率fs=100 Hz,h1=h2=h3=0.1,ω1、ω2、ω3分別為20,30,50 rad/s,噪聲強度D=1,Δω為0.2 rad/s。設定[ωmin,ωmax]為(0,60)。繪制待測信號的時域圖和頻域圖如圖5所示,波形顯示3 000個采樣點。 由于在實際工程問題中噪聲一般是不可忽略的,因此在時域圖和頻域圖中,可以輕易發(fā)現(xiàn)由于噪聲的污染,很難找到待測信號的頻率的存在。而相圖的檢測也同樣容易受到噪聲的干擾,如圖6所示。 將ω=20 rad/s的待測信號送入式(6)中的第一個方程,可以發(fā)現(xiàn)圖6(a)所顯示的相圖在噪聲的影響下,較難分辨是出于混沌狀態(tài)還是大周期態(tài)。而圖6(b)所顯示的則是在經(jīng)過尺度變換后消除了噪聲的影響,相圖是處于大周期態(tài)的。在下文中,都將以經(jīng)過尺度變換后的相圖闡述。 采用本文提出的優(yōu)化模型進行檢測,在ω=10 rad/s,20 rad/s,30 rad/s,40 rad/s,50 rad/s,60 rad/s,6個點處進行取值,如圖7所示。 (a)時域圖 (a)未經(jīng)過尺度變換的相平面軌跡 (b)經(jīng)過尺度變換的相平面軌跡 圖6 ω=20 rad/s的待測信號的相平面軌跡 通過實驗可以得知,圖7(b)、(c)、(e)所顯示的是分別為ω=20 rad/s、ω=30 rad/s和ω=50 rad/s的頻率分量所對應的相圖,相圖是明顯處于大周期態(tài)的。而圖7(a)、(d)、(f)所顯示的分別是在ω=10 rad/s、ω=40 rad/s和ω=60 rad/s的頻點所對應的相圖,相圖是處于混沌狀態(tài)的,符合實驗設定的預期條件。 實驗發(fā)現(xiàn): (1)未使用尺度變換,就很難將噪聲污染的相圖分辨出來,存在對相圖的誤判可能,尺度變換可以有效地消除噪聲對相圖的干擾,尺度系數(shù)R的存在有利于待測信號的提取。 (2)使用本方法成功檢測到了ω=20 rad/s,ω=30 rad/s,ω=50 rad/s的頻率分量。而在ω=60 rad/s這兩個頻點則相圖未發(fā)生狀態(tài)轉(zhuǎn)移,證明不存在該頻率的分量,具有普遍性和實用性。 (a)ω=10 rad/s的相平面軌跡 (b)ω=20 rad/s的相平面軌跡 (c)ω=30 rad/s的相平面軌跡 (d)ω=40 rad/s的相平面軌跡 (e)ω=50 rad/s的相平面軌跡 (f)ω=60 rad/s的相平面軌跡 (3)本方法實現(xiàn)的微弱信號檢測最低信噪比門限可以達到-24.2 dB,低于文獻[16]中提出的原檢測模型最低信噪比門限的-13.01 dB。 (4)Δω的選取如果過大,則會出現(xiàn)漏檢測現(xiàn)象;Δω的選取如果過小,則檢測效率就較低,因此需要找一個平衡點。一般采用比待測信號頻率低3個數(shù)量級的Δω可以實現(xiàn)限定頻率范圍內(nèi)全部頻率的檢測。本文所提出的方法在實際實驗中采用的Δω是0.2 rad/s,檢測到了3個待測信號的頻率。 (5)本文所提出的方法與Duffing振子并行處理方式的區(qū)別在于,本文所提出的方法是對頻域加帶通濾波器的基礎上進行而且僅用一個Duffing振子便可檢測到限定頻域內(nèi)全部的待測信號,而Duffing振子并行處理方式一般采用陣列的方式進行檢測,即需要較多振子。本文提出的方法節(jié)省了振子的數(shù)量,對運算條件需求較低,在實際實驗處理中也可以達到較為理想的效果。 本文研究了一種未知多頻微弱信號檢測方法,該方法針對原檢測模型只能同時檢測某一頻率的待測信號,克服了噪聲誘導對檢測結(jié)果的不利影響,提出了一種可同時檢測多個未知頻率微弱信號的Duffing振子優(yōu)化模型,可以實現(xiàn)提取未知多頻微弱信號各信號分量的頻率參數(shù),具有實用性和普遍性。仿真結(jié)果表明:該方法為微弱多頻信號處理提供了一種解決方案,拓寬了Duffing振子在微弱信號檢測領域的應用。2.2 待測信號提取模型
3 仿真與分析
4 結(jié)束語