孫偉瑋，陳俊喆，王 偉

一種基于短時傅里葉變換的羅蘭C信號頻譜分析方法

孫偉瑋1，陳俊喆2,3，王 偉2,3

（1 海軍裝備部；2 中國電子科技集團(tuán)公司第二十研究所，西安 710068；3 陜西省組合與智能導(dǎo)航重點(diǎn)實(shí)驗(yàn)室，西安 710068）

當(dāng)前對羅蘭C信號的頻譜分析基本上采用的是基于傅里葉變換法的頻譜分析方法，此類方法無法實(shí)現(xiàn)對信號頻譜的實(shí)時監(jiān)測。采用短時傅里葉變換法（STFT）對不同信干比、天波幅度和時延條件影響下的羅蘭C信號進(jìn)行理論建模與仿真，并建立信號的時頻數(shù)據(jù)庫。通過對 羅蘭C信號的STFT分析，并與時頻數(shù)據(jù)庫進(jìn)行譜峰比對，可有效解決實(shí)時性要求較高的天地波識別與天波時延測量問題。

羅蘭C；短時傅里葉變換；天波；頻譜分析

0 引言

羅蘭C是一種國際標(biāo)準(zhǔn)無線電定位導(dǎo)航和授時（Positioning，Navigation and Timing，PNT）服務(wù)系統(tǒng)，具有發(fā)射功率高、傳播距離遠(yuǎn)、相位穩(wěn)定性高等特點(diǎn)?？商峁?00 ns以上的授時服務(wù)和20 m的差分修正定位服務(wù)[1]。

隨著現(xiàn)代信號處理方式的引入，羅蘭C信號的處理方法已較多集中在頻域處理，如采用各種窗函數(shù)構(gòu)造的帶通濾波器是其抗噪聲的主要方式[2]；采用基于全相位快速傅里葉變換（all phase Fast Fourier Transform，apFFT）的陷波器針對點(diǎn)頻干擾進(jìn)行抑制[3]；采用基于快速傅里葉逆變換（Inverse Fast Fourier Transform，IFFT）的頻譜相除法[4]以及基于此法衍生的多信號分類（Multiple Signal Classification，MUSIC）算法[5]，用于解決天波干擾問題，對采樣數(shù)據(jù)進(jìn)行基于傅里葉變換法的頻域處理。由于采用一般的傅里葉變換法對采樣數(shù)據(jù)進(jìn)行處理的方式無法實(shí)時獲得與時間相關(guān)的頻域特征信息，同時處理存在滯后性，因此以往對羅蘭C信號的監(jiān)測較少進(jìn)行有關(guān)頻譜的實(shí)時監(jiān)測。

當(dāng)前，羅蘭C信號的天地波識別能力以及抗天波干擾能力已成為了衡量羅蘭C監(jiān)測接收機(jī)功能性能的一項(xiàng)重要指標(biāo)[6]。本文通過采用短時傅里葉變換法（Short Time Fourier Transform，STFT）對 羅蘭C信號進(jìn)行時頻分析，以應(yīng)對實(shí)時性要求較高的羅蘭C信號天地波識別和天波時延監(jiān)測。

1 羅蘭C信號頻域特性

1.1 羅蘭C信號頻域特性

式中，為與信號峰值幅度有關(guān)常數(shù)；為時間，單位μs；為包周差，單位μs；為系統(tǒng)相位編碼。

1.2 頻域波形

通過對羅蘭C信號進(jìn)行采樣獲得一段包含信號波形的采樣數(shù)據(jù)，對其進(jìn)行離散傅里葉變換（Discrete Fourier Transform，DFT），可以獲得其頻譜分布[8]，如式（2）所示：

當(dāng)采用1 MHz采樣率采集一段70 dBμV電平的典型羅蘭C信號數(shù)據(jù)，其頻譜分布如圖1藍(lán)色實(shí)線所示。

對包含有地波、天波及噪聲干擾的羅蘭C信號，可以表示為時域組合形式如式（3）所示：

典型的包含天波、地波和噪聲的羅蘭C信號頻譜分布如圖1紅色虛線所示。其中天波相對于地波延遲85 μs，幅度小2 dB，加30 dBμV的噪聲電平。

圖1 羅蘭C信號頻譜

從標(biāo)準(zhǔn)頻譜和加入天波的頻譜可以看出，即使不加入其他干擾，波形的頻譜也發(fā)生了變化，峰值不再在100 kHz，而形成了雙峰。通過對不同天波的時延、幅值等參數(shù)的變化研究，可以發(fā)現(xiàn)信號的頻譜圖也在發(fā)生各種變化，通過這種方法不能有效地總結(jié)出規(guī)律，很難用于識別羅蘭C信號及其天波的存在。

2 基于STFT的頻譜分析方法

基于傳統(tǒng)傅里葉變換獲得的信號頻譜是通過對一段采樣數(shù)據(jù)整體的頻域分布，從式（4）可以看出已經(jīng)丟失了時間關(guān)系，缺乏根據(jù)時間確定信號頻率的能力，對于許多需要進(jìn)行快速信號頻譜分析的場合難以使用。

3）可以看作用基函數(shù)，如式（6）所示：

采用STFT方法對滿足式（1）的一段標(biāo)準(zhǔn)信號和滿足式（3）的附加天波和噪聲的信號進(jìn)行處理，其頻譜圖圖2和圖3所示，采用的是5 μs的軸移動步進(jìn)和30 μs的分段進(jìn)行處理。從圖中可以直觀地看出天波在地波后多長時間產(chǎn)生，并如何發(fā)生的變化。

圖2 標(biāo)準(zhǔn)信號STFT頻譜

圖3 加天波和噪聲的信號STFT頻譜

在此基礎(chǔ)上，可對每一次STFT處理后的頻譜數(shù)據(jù)在以中心頻率為100kHz的附近搜索頻率峰值點(diǎn)，如圖4所示。仿真設(shè)置的地波前沿起始時刻是150 μs處，實(shí)測峰值出現(xiàn)在150 μs處；仿真設(shè)置的天波起始點(diǎn)在地波后85 μs（即235 μs處）出現(xiàn)，實(shí)測天波在地波后110 μs處出現(xiàn)（即260 μs處）出現(xiàn)明顯變化，較實(shí)際天波出現(xiàn)位置滯后約25 μs。

圖4 頻譜峰值位置變化規(guī)律

可以看出，STFT方法對波形前沿的檢測十分靈敏，一般噪聲譜峰值點(diǎn)的分布為隨機(jī)分布，而信號譜的峰值分布主要集中在100 kHz附近，信號譜的功率較高，容易識別出最高峰值；同時天波在出現(xiàn)時刻前后，信號譜峰值中心頻率會發(fā)生波動，波動的起始時刻一般在天波的起始點(diǎn)附近。這對于研究天波信號的起始點(diǎn)位置具有重要意義。

3 仿真驗(yàn)證

3.1 信噪比影響測試

假設(shè)不存在干擾信號，僅存在高斯白噪聲信號，當(dāng)改變輸入信號的信噪比時，通過采用STFT法對其進(jìn)行頻域信號識別，通過搜索中心頻率附近連續(xù)峰值頻率的起始位置，獲得信號被識別的起始時間，包括地波起始時間和天波起始時間。測試采用5 μs移動步進(jìn)。測試結(jié)果如表1所示。

表1 不同信噪比條件下起始時刻

通過多次測試可以看出，在正信噪比條件下，STFT法可以識別信號脈沖的前沿起始位置，但是隨著信噪比的下降，識別精度越來越差，對于天波下降更為明顯。因此在數(shù)據(jù)進(jìn)行STFT處理前需要進(jìn)行有效的帶通濾波，這樣會對識別的準(zhǔn)確性有很大提升。

3.2 天波幅值影響測試

假設(shè)在相同的信噪比和天波時延的情況下，改變天波的幅值，測試STFT法對信號識別的影響。設(shè)置地波信號70 dBμV，天波時延85 μs，信噪比為20 dB的測試條件，天波的幅值從40 ～100 dBμV區(qū)間，地波與天波的信干比（Signal to Interference Ratio，SIR）為30 dB～-30 dB，10 dB步進(jìn)。

這里首先對包含天波干擾的不同的信號波形進(jìn)行頻譜分析，如圖5所示。可以看出不同的條件下信號的頻譜是有很大差異的，當(dāng)>10 dB時，信號譜分量主要以地波為主；當(dāng)在0 dB附近時，信號譜存在明顯的中心頻率偏移或雙峰譜的現(xiàn)象；當(dāng)<-10 dB時，信號譜分量主要以天波為主。由于實(shí)際信號地波場強(qiáng)與天波場強(qiáng)的關(guān)系無法直接得知，沒有辦法直接從一般傅里葉變換中分析出地波頻譜和天波頻譜。從圖5這種直接傅里葉變換的頻譜分析方法中很難直觀地反映出信號頻譜隨時間的變化特性。為了更好地觀察包含天波的信號頻譜特性，采用STFT法進(jìn)行更詳細(xì)的分析。

圖5 不同信干比下的天地波信號譜

幾種不同條件下的STFT法的天地波信號頻譜分布如圖6所示。其中橫軸為頻率，縱軸為分段所在的時間。不同于圖2和圖3的3D繪圖表示法，圖6采用2D平面表示，用不同的顏色深度來表示信號的功率幅值，越紅則幅值越高，越藍(lán)幅值越低?？梢钥闯鲭S著天波信號幅度的增強(qiáng)，天波信號在時間軸上的頻譜分量顯示得越清晰；天波信號幅度與地波信號幅度越接近，且天地波的分界越明顯，越容易識別。

3.3 天波時延影響測試

以往的研究中，對于天波時延的估計(jì)主要針對在地波起始點(diǎn)后42.5～160 μs時間段內(nèi)出現(xiàn)的一跳天波信號。通過對在這一時延范圍內(nèi)的全部頻譜進(jìn)行仿真，可以建立一套基于時延的頻譜數(shù)據(jù)庫。時延45 μs選圖如圖7所示。

圖7 時延45 μs下STFT法的天地波信號譜

為了與下一節(jié)實(shí)測數(shù)據(jù)進(jìn)行比對。通過對實(shí)際信號接收的數(shù)據(jù)與數(shù)據(jù)庫的數(shù)據(jù)進(jìn)行比對，采用對譜峰相對位置進(jìn)行測量比較的方法，可以在頻域識別該接收信號的地波起始點(diǎn)和天波時延等信息。

通過采用對模擬仿真信號進(jìn)行在不同信噪比條件下的地波起始時間、天波時延時間測試，并對不同天波幅值和天波時延條件下的頻譜進(jìn)行STFT仿真，并建立比對頻譜數(shù)據(jù)庫。該方法對信號處理的實(shí)時性要求較高。通過對仿真計(jì)算統(tǒng)計(jì)，該處理方法可在μs量級處理完成，完全可以用于實(shí)時信號的快速識別與分析。

4 實(shí)測驗(yàn)證

為了驗(yàn)證該方法對羅蘭C信號的識別有效性，通過在西安地區(qū)接收宣城臺（8390 M）的羅蘭C信號，采用STFT法對接收信號進(jìn)行分析，獲得信號起始時間和天波時延量，反算電離層高度，并與理論高度進(jìn)行比對。測試本地時間為6月8日17時20分左右，該時間段為日落之前時段，電離層較低但活動較為劇烈。測試結(jié)果如圖8和圖9所示。

圖9 宣城臺天地波信號譜峰值

通過圖8對信號譜與數(shù)據(jù)庫的信號譜進(jìn)行比對，其天地波譜峰分離程度與時延45 μs的信號譜較為接近，通過對圖9譜峰最大值位置的分析，信號的起始時間約在120～130 μs之間，譜峰最大偏移量在200 μs處出現(xiàn)，根據(jù)之前的仿真結(jié)論，天波出現(xiàn)時刻約在175 μs處，因此天波的相對地波的時延量約為45～55 μs。綜合上述判斷天波時延應(yīng)該在45 μs處附近。根據(jù)西安到宣城臺的大地球面距離約 1 000 km，根據(jù)長波天波傳播路徑的分析和研究[9]的結(jié)論，可以計(jì)算出電離層的高度約在60～65 km，與D電離層在測試時的本地時間下中緯度地區(qū)的電離層高度基本一致，STFT的結(jié)果基本得到了印證。

5 結(jié)語

本文通過分析羅蘭C的頻域特性，提出了一種基于STFT進(jìn)行羅蘭C信號頻譜分析的方法。該方法具有實(shí)時性高、信號譜特征識別能力強(qiáng)、對不同信噪比和天波時延信號由較高的分辨能力等優(yōu)點(diǎn)。通過對該方法進(jìn)行的仿真和試驗(yàn)驗(yàn)證，證明了這些優(yōu)點(diǎn)，可用于未來羅蘭C信號的監(jiān)測、天波時延監(jiān)測，甚至電離層高度監(jiān)測，若羅蘭C用戶接收機(jī)采用此方法，可以實(shí)現(xiàn)天地波快速識別；而且，在地波信號覆蓋區(qū)外，可利用天波信號實(shí)現(xiàn)遠(yuǎn)距離授時。

[1] WENHE YAN，KUNJUAN ZHAO，SHIFENG LI，et al. Precise Loran-C Signal Acquisition Based on Envelope Delay Correlation Method[J]. MDPI Sensors，2020，2329（20）：1-15.

[2] D. LAST，Y. BIAN. Carrier wave interference and Loran-C receiver performance[J]. IEEE PROCEEDINGS——F，1993，140（5）：273-283.

[3] 林洪文，周旻，朱四華，等. 基于apFFT的羅蘭C窄帶干擾頻域抑制方法[J]. 海軍航空工程學(xué)院學(xué)報(bào)，2012，27（6）：613-617.

[4] A. MOHAMMED，D. LAST. IFFT technique for skywave detection in Loran-C receivers[J]. Electronic Letters，2001，37（6）：398-400.

[5] Y. BIAN，J. D. Last. Loran-C skywave delay estimation using eigen-decomposition techniques[J]. Electronic Letters，1995，31（2）：133-134.

[6] 潘峰，李國俊，楊大峰. 長河二號授時監(jiān)測接收機(jī)的設(shè)計(jì)與實(shí)現(xiàn)[J]. 宇航計(jì)測技術(shù)，2020，40（3）：1-5.

[7] HONGLEI QIN，XIAOQIN JIN，CONG LI，et al. MEDLL-based method of ground-wave and cycle identification for Loran-C signal[J]. 14th IEEE International Conference on Electronic Measurement & Instruments，2019：114-123.

[8] 胡廣書. 數(shù)字信號處理——理論、算法與實(shí)踐[M]. 北京：清華大學(xué)出版社，2003.

[9] 王述香，陳秀榮，辛永訓(xùn)，等. 天波的多跳傳播模型研究[J]. 數(shù)學(xué)建模及其應(yīng)用，2018，7（3）：49-57.

Method for Spectrum Analysis of Loran-C Signals Based on Short Time Fourier Transform

SUN Weiwei, CHEN Junzhe, WANG Wei

The current spectrum analysis method of Loran-C signals based on Fourier transform method, which can not realize the real-time monitoring of the signal spectrum. The Short-Time Fourier Transform (STFT) method is used to model and simulate Loran-C signals under different signal-to-interference ratio, sky wave amplitude and delay conditions, and the time-frequency database of the signal is established. The STFT analysis of Loran-C signal and spectral peak comparison with time-frequency database can effectively solve the problem of sky wave identification and sky wave delay measurement with high real-time requirements.

Loran-C; Short-Time Fourier Transform; Sky Wave; Spectral Analysis

TN961

A

1674-7976-(2022)-05-339-05

2022-08-01。孫偉瑋（1985.12—），山東平度人，碩士研究生，工程師，主要研究方向?yàn)閷?dǎo)航專業(yè)裝備質(zhì)量監(jiān)督。