程凌峰,倪淑燕,陳世淼,廖育榮

(航天工程大學(xué) a.研究生院;b.電子與光學(xué)工程系,北京 101416)

0 引 言

近幾年來,Starlink、OneWeb等低軌互聯(lián)網(wǎng)星座快速發(fā)展,已經(jīng)成為一種新型通信方式,具有廣泛的應(yīng)用前景。對(duì)于飛機(jī)、導(dǎo)彈等用戶,由于低軌衛(wèi)星高速運(yùn)動(dòng),終端也在高動(dòng)態(tài)運(yùn)動(dòng),在這種情況下接收信號(hào)受到較強(qiáng)多普勒效應(yīng)的影響[1],體現(xiàn)為較大的多普勒頻偏及高階頻偏變化率[2]。由于衛(wèi)星軌道可預(yù)報(bào),因此由衛(wèi)星運(yùn)動(dòng)引起的多普勒頻移可進(jìn)行預(yù)補(bǔ)償。但終端動(dòng)態(tài)性引起的多普勒頻移則需要利用接收機(jī)進(jìn)行捕獲,當(dāng)接收信號(hào)時(shí)間內(nèi)頻偏變化超出快速傅里葉變換(Fast Fourier Transform,FFT)分辨率的情況下,傳統(tǒng)捕獲方法[3-4]難以對(duì)信號(hào)進(jìn)行有效的估計(jì)。

通信系統(tǒng)在設(shè)計(jì)幀結(jié)構(gòu)的時(shí)候,通常在每一幀的最前端添加導(dǎo)頻序列用于捕獲,該序列在經(jīng)過高動(dòng)態(tài)信道之后變?yōu)檎{(diào)頻信號(hào)。目前針對(duì)調(diào)頻信號(hào)的參數(shù)估計(jì)方法主要有兩類：參數(shù)域類方法和時(shí)頻分析法。參數(shù)域類方法比較常用的有調(diào)頻率逼近法[5-6]和降階類參數(shù)估計(jì)算法[7-8],調(diào)頻率逼近法對(duì)線性調(diào)頻信號(hào)有較好的估計(jì)性能,但難以兼顧精度與復(fù)雜度,且不適用于高階調(diào)頻信號(hào);降階類參數(shù)估計(jì)算法可用于非線性調(diào)頻信號(hào),但是存在高階參數(shù)估計(jì)量對(duì)低階參數(shù)估計(jì)時(shí)誤差積累的問題,同時(shí)可能存在偽峰,抗噪性較差。時(shí)頻分析法是分析非平穩(wěn)信號(hào)的理想工具[9],它對(duì)信號(hào)進(jìn)行時(shí)域和頻域聯(lián)合分析獲得時(shí)頻圖,通過提取時(shí)頻圖脊線即瞬時(shí)頻率(Instantaneous Frequency,IF)實(shí)現(xiàn)對(duì)頻率的粗估,如文獻(xiàn)[10]就利用分?jǐn)?shù)階傅里葉變換實(shí)現(xiàn)了對(duì)線性調(diào)頻信號(hào)的頻率估計(jì)。對(duì)于高階調(diào)頻信號(hào)的估計(jì),過去通常使用短時(shí)傅里葉變換和小波變換,但是受到海森堡(Heisenberg)測(cè)不準(zhǔn)原理的限制[11],其時(shí)頻分布是模糊的,不利于脊線的提取。為了提高傳統(tǒng)方法的能力,近幾十年來學(xué)者們提出了時(shí)頻重排(Reassignment Method,RM)[12]、同步壓縮(Synchrosqueezing Transform,SST)[13]、高階同步壓縮(High-Order Synchro Squeezing Transform,HOSST)[14]和同步提取變換(Synchro Extracting Transform,SET)[15]等方法。時(shí)頻重排法具有較好的時(shí)頻聚集能力,但是無法實(shí)現(xiàn)信號(hào)重構(gòu);同步壓縮法對(duì)快時(shí)變信號(hào)的聚集能力較差,存在頻率估計(jì)誤差逐步增大的問題;高階同步壓縮對(duì)噪聲魯棒性較差,低信噪比下時(shí)頻圖發(fā)散嚴(yán)重;同步提取變換方法是短時(shí)傅里葉變換的后處理形式,去除大部分時(shí)頻點(diǎn),提高時(shí)頻能量聚集度,在重構(gòu)信號(hào)時(shí)減少了計(jì)算量,但是存在對(duì)快時(shí)變信號(hào)魯棒性較差的問題[16]。

針對(duì)以上問題,本文根據(jù)高動(dòng)態(tài)信號(hào)特點(diǎn)將其建模為非線性調(diào)頻信號(hào),并提出一種基于多重同步壓縮(Multisynchro Squeezing Transform,MSST)的時(shí)頻分析法提高時(shí)頻聚焦能力。在此基礎(chǔ)上對(duì)傳統(tǒng)脊線提取方法進(jìn)行了改進(jìn),設(shè)計(jì)了一種分段能量泛函最小化法,解決了由于噪聲導(dǎo)致的提取點(diǎn)錯(cuò)誤問題,提高低信噪比下信號(hào)估計(jì)能力。

1 高動(dòng)態(tài)信號(hào)建模

發(fā)射信號(hào)經(jīng)過信道傳輸進(jìn)入接收機(jī),并經(jīng)過自動(dòng)增益控制、下變頻和匹配濾波后,進(jìn)入載波同步。信號(hào)模型可表示為

(1)

式中：s(k)表示調(diào)制數(shù)據(jù),k為整數(shù);gT(t)為信號(hào)成形濾波器的沖激響應(yīng);T為符號(hào)周期;Δf和Δθ為接收機(jī)接收到的信號(hào)載波頻率、載波相位與本地載波頻率、相位的差值;w(nT)是均值為0的高斯白噪聲。

導(dǎo)頻序列通常為無數(shù)據(jù)變化的常數(shù)序列：

s(k)=1,k=0,1,…,P-1

(2)

式中：P為導(dǎo)頻長。

暫不考慮成形濾波器的拖尾,則經(jīng)過下變頻和濾波后的基帶信號(hào)可表示為

(3)

受相對(duì)速度和相對(duì)速度變化率等的影響,可將時(shí)變基帶信號(hào)相位近似高階泰勒級(jí)數(shù)展開：

(4)

式中：f0,f1,f2分別為頻偏、一階頻偏變化率和二階頻偏變化率。再考慮到導(dǎo)頻數(shù)據(jù)較短的情況,忽略了更高階頻偏變化率分量造成的影響,從而載波捕獲的問題轉(zhuǎn)換為非線性調(diào)頻信號(hào)的參數(shù)估計(jì)問題。

2 基于多重同步壓縮的高動(dòng)態(tài)信號(hào)時(shí)頻分析

2.1 同步壓縮基本原理

同步壓縮是基于短時(shí)傅里葉變換(Short Time Fourier Transform,STFT)進(jìn)行后處理得到的。對(duì)于信號(hào)r(t),STFT可表示為

(5)

式中：g(t)為窗函數(shù);r(t)為輸入信號(hào)。受到海森堡測(cè)不準(zhǔn)原理的制約,STFT無法同時(shí)實(shí)現(xiàn)時(shí)頻的高分辨率。SST方法實(shí)際上是對(duì)STFT在頻率方向進(jìn)行能量重新分配的一種方法,如圖1所示,主要用于純諧波信號(hào)分析,其過程可以表示為[17]

圖1 同步壓縮方法時(shí)頻能量分配Fig.1 Time-frequency energy distribution using synchronous compression method

(6)

接收信號(hào)的SST表達(dá)式為

(7)

(8)

(9)

2.2 高動(dòng)態(tài)信號(hào)的同步壓縮

對(duì)于本文研究的高動(dòng)態(tài)信號(hào),由于STFT的窗函數(shù)較短,二階頻偏變化率影響較小,因此可以將其所加窗內(nèi)信號(hào)近似為線性調(diào)頻信號(hào),即

r(τ)=Aej[φ(t)+2πf0(τ-t)+πf1(τ-t)2]

(10)

同時(shí),定義窗函數(shù)為高斯窗g(t)=e-0.5t2,則

(11)

(12)

因此將式(12)兩邊對(duì)時(shí)間求偏導(dǎo)可得

(13)

因?yàn)槭?13)是復(fù)數(shù)形式,不能直接用于計(jì)算,取實(shí)部可得

(14)

2.3 高動(dòng)態(tài)信號(hào)的多重同步壓縮

為了提高處理快時(shí)變信號(hào)時(shí)的頻率分辨率,前人通過采用高階的同步壓縮方法以應(yīng)對(duì)多階的瞬時(shí)頻率,但是高階同步壓縮在較低信噪比環(huán)境下性能惡化嚴(yán)重。為此,本文設(shè)計(jì)了一種高動(dòng)態(tài)信號(hào)的多重同步壓縮方法。

(15)

進(jìn)一步可得

(16)

圖2 MSST頻率壓縮過程Fig.2 MSST frequency compression process

(17)

(18)

經(jīng)過多次迭代,模糊的STFT時(shí)頻分布能量得到集中,瞬時(shí)頻率的估計(jì)會(huì)越來越接近信號(hào)的真實(shí)頻率,在脊線提取環(huán)節(jié)就會(huì)獲得更為準(zhǔn)確的瞬時(shí)頻率脊線,提高參數(shù)估計(jì)精度。

3 基于分段能量泛函最小化的瞬時(shí)頻率估計(jì)

時(shí)頻分析為信號(hào)參數(shù)估計(jì)提供了有力的工具,瞬時(shí)頻率作為信號(hào)的重要物理量,其從時(shí)頻圖提取的準(zhǔn)確性對(duì)非平穩(wěn)信號(hào)的研究具有重要的意義?；跁r(shí)頻分析的瞬時(shí)頻率提取方法有很多,如譜峰檢測(cè)法[18]、能量泛函最小化法[19]等。譜峰檢測(cè)法在高信噪比環(huán)境下性能較好,但是低信噪比下噪聲會(huì)導(dǎo)致時(shí)頻圖譜峰點(diǎn)的偏移,造成較大誤差。

3.1 能量泛函最小化法

能量泛函最小化法對(duì)信號(hào)時(shí)頻分析的結(jié)果增加懲罰函數(shù),選取匹配參數(shù),計(jì)算出使得曲線平滑且能量最大的瞬時(shí)頻率曲線。其計(jì)算公式為

(19)

式中：λ和β是調(diào)整正則化水平的兩個(gè)參數(shù)。通過式(19)可以看出,時(shí)頻因子的能量以及提取瞬時(shí)頻率的平滑度都是要被考慮的因素。文獻(xiàn)[20]中提出的瞬時(shí)頻率提取算法的起始點(diǎn)選擇為全域最大值,這種方法的優(yōu)點(diǎn)是可以明顯減小計(jì)算量,但是在某些情況下不能保證全局最大值精確地位于真實(shí)的瞬時(shí)頻率上,這就偏離了式(19)的初衷,如異常情況引起的較大波動(dòng)或者信噪比較低時(shí)導(dǎo)致的不正確起始點(diǎn)。

3.2 分段能量泛函最小化法

為了提高瞬時(shí)頻率提取過程的效率和可靠性,考慮將整個(gè)時(shí)頻平面分割成幾個(gè)區(qū)段。具體步驟如下：將給定的時(shí)頻平面Ts(h,m)分為F個(gè)部分,設(shè)置最大容許頻率變化為Δ,數(shù)據(jù)長度為L。對(duì)劃分的時(shí)頻平面依次分段確定起始點(diǎn),并進(jìn)行前向后向同時(shí)搜索瞬時(shí)頻率脊線,最后根據(jù)時(shí)頻能量從這些軌跡中選出最佳的瞬時(shí)頻率軌跡。

結(jié)合時(shí)頻分析和頻率提取過程,高動(dòng)態(tài)載波頻率捕獲算法的整體流程如圖3所示。

圖3 高動(dòng)態(tài)捕獲方法流程Fig.3 Flowchart of high dynamic aquisition method

4 性能仿真分析

4.1 仿真條件

由于衛(wèi)星運(yùn)動(dòng)軌道確定,因此由衛(wèi)星運(yùn)動(dòng)造成的多普勒頻移一般會(huì)采用預(yù)補(bǔ)償方式,高動(dòng)態(tài)性主要考慮由終端引起的。美國JPL提出了一種具有典型性的高動(dòng)態(tài)模型,假設(shè)飛行器運(yùn)動(dòng)速度最大為2 km/s(6Ma),加速度最大為400 m/s2(40g),加加速度最大為1 000 m/s3(100g/s)。將接收信號(hào)建模為二階非線性調(diào)頻信號(hào),按移動(dòng)衛(wèi)星通信常采用的S頻段(3 GHz)計(jì)算,多普勒頻移范圍±20 kHz,變化率范圍±4 kHz/s,二階變化率范圍±10 kHz/s2。采樣頻率設(shè)為fs=50 kHz,數(shù)據(jù)長度N=5 000。頻偏、一階調(diào)頻率、二階調(diào)頻率在取值范圍內(nèi)隨機(jī)取值,在計(jì)算頻率估計(jì)均方根誤差和捕獲概率時(shí)采用500次蒙特卡羅仿真,比較STFT、RM、SST、SST3、SET和MSST這幾種時(shí)頻分析方法在高動(dòng)態(tài)非線性調(diào)頻信號(hào)參數(shù)估計(jì)中的性能。其中,SST3表示三階同步壓縮方法(Third-order Synchro Squeezing Transform),是高階同步壓縮的典型方式。

4.2 迭代次數(shù)的確定

作為一種迭代算法,需要確定MSST迭代終止的條件。本文采用時(shí)頻分析瑞利熵和重構(gòu)信號(hào)的信噪比變化來確定合適的迭代次數(shù)。

時(shí)頻聚集性是時(shí)頻分析方法性能優(yōu)劣的重要指標(biāo)之一,通過將時(shí)頻分析和信息熵結(jié)合,提出利用瑞利熵來表征時(shí)頻分析方法性能。α階瑞利熵的定義為[21]

(20)

Ryα越小,時(shí)頻聚集程度越高。瑞利熵可以確保MSST時(shí)頻聚集到合適的水平,當(dāng)?shù)^程中連續(xù)時(shí)頻分析結(jié)果的瑞利熵沒有明顯的變化時(shí),迭代過程可以中止。

對(duì)時(shí)頻分析結(jié)果進(jìn)行重構(gòu),比較不同迭代次數(shù)重構(gòu)后信號(hào)的信噪比,如果相鄰兩次信噪比不再明顯變化,迭代過程可以中止。圖4和圖5分別為迭代次數(shù)對(duì)瑞利熵和重構(gòu)的影響,可以看出,迭代次數(shù)越多,瑞利熵越小,重構(gòu)后的信號(hào)信噪比越高,經(jīng)過3次迭代,瑞利熵和信噪比趨于穩(wěn)定。除了MSST的重構(gòu)信噪比外,圖5中還給出了STFT、SST、SST3和SET的重構(gòu)信噪比,輸入信噪比為5 dB。可以看出,除MSST外其他幾種后處理方法都會(huì)實(shí)現(xiàn)一定的信噪比提升,但是STFT重構(gòu)時(shí)會(huì)降低信噪比。SST略有提升(MSST一次迭代時(shí)),SST3大約提升4 dB,SET提升約6 dB,而MSST在穩(wěn)定后信噪比提升了約9 dB,性能提升非常明顯,這證明了后處理技術(shù)的優(yōu)勢(shì)。MSST的這種能力可以考慮用于信號(hào)的信噪比增強(qiáng)。

圖4 迭代次數(shù)對(duì)瑞利熵的影響Fig.4 Influence of iteration times on Rayleigh entropy

圖5 迭代次數(shù)對(duì)重構(gòu)信噪比的影響Fig.5 Influence of iteration times on the reconstruction signal-to-noise ratio

4.3 時(shí)頻聚集能力

假設(shè)信號(hào)信噪比為5 dB,STFT、RM、SST、SST3、SET和MSST這6種方法的時(shí)頻二維能量分布如圖6所示,其中右側(cè)為時(shí)頻分析局部放大圖。對(duì)除了STFT之外的時(shí)頻圖采用統(tǒng)一刻度進(jìn)行比較可以看出,STFT的能量聚集最差,局部放大后幾乎難以分辨脊線;SST的瞬時(shí)頻率受到估計(jì)模型和噪聲的限制,時(shí)頻發(fā)散,時(shí)頻點(diǎn)能量較低;RM方法由于是從時(shí)間和頻率兩個(gè)方向重新分配時(shí)頻能量,因此得到了比SST更集中的時(shí)頻圖,但放大后時(shí)頻圖的脊線還是處于發(fā)散狀態(tài);SST3雖然是針對(duì)二階非線性調(diào)頻信號(hào)設(shè)計(jì),相對(duì)于SST方法模型更為準(zhǔn)確,但是在較低信噪比下性能欠佳,時(shí)頻圖較為發(fā)散;SET方法通過去除模糊的時(shí)頻點(diǎn)實(shí)現(xiàn)能量集聚,在某些時(shí)段聚集性能優(yōu)于MSST,但由于對(duì)快時(shí)變的魯棒性欠佳,在時(shí)頻線的交界處能量斷層且在噪聲的影響下時(shí)頻線出現(xiàn)上下波動(dòng)的現(xiàn)象,造成估計(jì)誤差增大;MSST實(shí)現(xiàn)了最集中的時(shí)頻分析結(jié)果,能夠看到清楚的時(shí)頻圖脊線,為后續(xù)提取脊線進(jìn)行準(zhǔn)確瞬時(shí)頻率估計(jì)提供了保證。

圖6 不同時(shí)頻分析方法分析結(jié)果及局部放大圖Fig.6 Analysis results and local magnified images using different time-frequency analysis methods

用瑞利熵對(duì)6種方法的性能進(jìn)行定量分析,結(jié)果如圖7所示?？梢钥闯?輸入信噪比在-10～15 dB時(shí),MSST算法的瑞利熵遠(yuǎn)小于其他幾種方法,這也證明了在幾種方法中MSST的時(shí)頻聚焦能力最優(yōu)。

圖7 不同時(shí)頻分析方法的瑞利熵隨信噪比變化情況Fig.7 Variation of Rayleigh entropy with signal-to-noise ratio using different time-frequency analysis methods

4.4 頻率估計(jì)精度

圖8為不同方法的頻率估計(jì)均方根誤差隨信噪比變化情況,可見MSST方法的頻率估計(jì)均方根誤差始終是最小的,在較高信噪比(大于10 dB)時(shí),各個(gè)方法的均方根誤差變化趨于平穩(wěn),MSST方法的均方根誤差約為STFT方法的1/3,且均方根誤差大小與時(shí)頻聚集能力成負(fù)相關(guān)。在較低信噪比(小于0 dB)時(shí),雖然MSST的均方根誤差依然是最小的,但是其他5種方法的均方根誤差不是完全按照時(shí)頻聚集能力的順序排列的,這說明作為STFT的后處理方法,不僅將有用信號(hào)的瞬時(shí)頻率時(shí)頻能量聚集,同時(shí)也會(huì)將噪聲聚集,在瞬時(shí)頻率估計(jì)的時(shí)候造成過大的誤差。很明顯,MSST對(duì)信號(hào)的時(shí)頻聚集能力強(qiáng)于噪聲聚集能力,而SST、RM和SST3在低信噪比時(shí)瞬時(shí)頻率的時(shí)頻能量嚴(yán)重發(fā)散,無法實(shí)現(xiàn)準(zhǔn)確地瞬時(shí)頻率估計(jì),甚至相較于STFT的估計(jì)精度都出現(xiàn)下滑,其中SST3的頻率估計(jì)誤差在信噪比-1 dB以下保持在高位,噪聲魯棒性較差,但是隨著信噪比的提升,SST3能夠準(zhǔn)確分析快時(shí)變信號(hào)的能力逐漸突顯,在較高信噪比下頻率估計(jì)誤差僅次于MSST。SET方法在較低信噪比時(shí)頻率估計(jì)性能較好,但是在高信噪比時(shí)僅略優(yōu)于STFT,性能整體表現(xiàn)平穩(wěn)。MSST的這個(gè)特性對(duì)于在較低信噪比下的載波捕獲有重要的意義。本文中分段能量泛函最小化法的分段數(shù)取為400。由圖8可知,信噪比在-5 dB以下時(shí)分段意義重大,即使是MSST方法,在這種很低的信噪比下如果不分段頻率估計(jì)性能下滑都極其嚴(yán)重,在-10 dB時(shí)頻率估計(jì)均方根誤差增大了近10倍,在-5 dB以上分段與否對(duì)性能影響很小,可以根據(jù)實(shí)際應(yīng)用場(chǎng)景選擇合理的分段數(shù)實(shí)現(xiàn)精度和運(yùn)算量的平衡。

4.5 捕獲概率與時(shí)間

圖9比較了不同算法的捕獲概率。捕獲概率通常定義為信號(hào)頻率估計(jì)誤差小于20 Hz的概率。MSST方法在捕獲概率收斂后可達(dá)到99%,且在低信噪比條件下仍能保持95%以上,明顯優(yōu)于其他方法。同時(shí)也能發(fā)現(xiàn)SST、RM、SST3和SET這4種方法在高信噪比條件下捕獲概率均收斂到95%以上,而STFT的捕獲概率收斂到80%左右,也證明了對(duì)STFT后處理的必要性。同樣也可以發(fā)現(xiàn),在-6 dB以下頻率估計(jì)時(shí)不分段的缺點(diǎn)逐漸明顯,在-10 dB時(shí)MSST法的捕獲概率降到90%以下。

圖9 不同方法捕獲概率隨信噪比變化情況Fig.9 Variation of aequisition probability with signal-to-noise ratio using different methods

在實(shí)際應(yīng)用時(shí),時(shí)頻分析方法的運(yùn)行效率非常關(guān)鍵,這決定了其能否應(yīng)用于實(shí)際工程中。假設(shè)時(shí)間點(diǎn)數(shù)和頻率點(diǎn)數(shù)分別為Nt和Nf,STFT的時(shí)間計(jì)算復(fù)雜度為O(NtNflb(Nf),后處理方法的多重壓縮、求導(dǎo)及重新分配時(shí)頻能量等步驟的計(jì)算復(fù)雜度為低階項(xiàng),可忽略。采用測(cè)試的計(jì)算機(jī)配置為i7-10750H 2.6 GHz,Win10系統(tǒng),Matlab版本為R2020a。所需時(shí)間如表1所示,所有方法都在15 s內(nèi)完成處理。MSST相對(duì)于SST方法時(shí)間增加較少,且隨著迭代次數(shù)的大幅增加,MSST方法所需計(jì)算時(shí)間增加幅度不大,即使選擇了較大的迭代次數(shù),MSST方法的計(jì)算時(shí)間仍明顯小于SST3,具有較好的實(shí)際使用價(jià)值。

表1 不同時(shí)頻分析方法所需計(jì)算時(shí)間Tab.1 Calculation time required for different time-frequency analysis methods

5 結(jié)束語

本文針對(duì)傳統(tǒng)的時(shí)頻圖頻率估計(jì)方法起始點(diǎn)選取錯(cuò)誤導(dǎo)致脊線提取異常的問題,提出在時(shí)頻平面上分段選取起始點(diǎn),前向后向同時(shí)搜索,獨(dú)立提取脊線,提高頻率估計(jì)準(zhǔn)確度。在高信噪比條件下頻率估計(jì)均方根誤差趨于收斂時(shí),MSST方法的均方根誤差約為STFT和SET的1/3,為SST、RM和SST3等后處理方法的1/2,捕獲概率為99%以上,且在低信噪比條件下,MSST的頻率估計(jì)性能波動(dòng)不大,捕獲概率為95%以上,取得了較好捕獲效果。同時(shí),在確定MSST最佳迭代次數(shù)時(shí)發(fā)現(xiàn),MSST存在提升重構(gòu)信號(hào)信噪比的能力,也為增強(qiáng)信號(hào)提供了一種思路,提高了MSST的實(shí)用性。

MSST方法是在SST方法框架下的方法,SST方法是假設(shè)所取窗內(nèi)信號(hào)為單頻諧波信號(hào),雖然通過多步壓縮提高時(shí)頻聚集度,但在本文場(chǎng)景下還是會(huì)導(dǎo)致頻率估計(jì)誤差隨時(shí)間積累。SST3是一種很好地技術(shù),它就是針對(duì)信號(hào)是二階非線性調(diào)頻信號(hào)的基礎(chǔ)上,但是受制于本身較差的噪聲魯棒性,在較低信噪比條件下性能很差,后續(xù)可以嘗試對(duì)其改進(jìn)用于信號(hào)捕獲。