含噪多頻信號(hào)頻率估計(jì)算法

陳 鵬， 劉春華， 蘇 欣， 涂亞慶， 趙少美

(1. 中國(guó)空氣動(dòng)力研究與發(fā)展中心 空氣動(dòng)力學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室，四川 綿陽(yáng) 621000；2. 中國(guó)空氣動(dòng)力研究與發(fā)展中心 設(shè)備設(shè)計(jì)及測(cè)試技術(shù)研究所，四川 綿陽(yáng) 621000；3. 陸軍勤務(wù)學(xué)院 軍事物流系，重慶 401311)

多頻信號(hào)可理解為多個(gè)單頻信號(hào)的疊加，其頻率估計(jì)廣泛應(yīng)用于線性系統(tǒng)識(shí)別、低頻機(jī)械光譜學(xué)、電力系統(tǒng)、核磁共振波譜分析以及無(wú)損檢測(cè)等諸多領(lǐng)域[1-3]。

目前，針對(duì)多頻信號(hào)的頻率估計(jì)算法大致可分為兩大類：時(shí)域法和頻域法。時(shí)域法如線性預(yù)測(cè)法(linear prediction, LP)[4]、計(jì)算量更加復(fù)雜的多重信號(hào)分類法(multiple signal classification, MUSIC)[5]、旋轉(zhuǎn)不變估計(jì)法(estimating signal parameter via rotational invariance techniques, ESPRIT)[6]。時(shí)域法分辨率高，但所需采樣序列長(zhǎng)，導(dǎo)致計(jì)算量大，且抗噪性差，不利于實(shí)際應(yīng)用。頻域法主要利用離散傅里葉變換法(discrete Fourier transform，DFT)對(duì)信號(hào)進(jìn)行頻譜分析，可利用DSP(digital signal processing)和FPGA(field programmable gate array)等硬件直接實(shí)現(xiàn)、計(jì)算速度快，抗噪性好，是研究多頻信號(hào)頻率估計(jì)的主要方向[7]。

單頻信號(hào)是多頻信號(hào)的一種特例，也是研究重點(diǎn)。文獻(xiàn)[8]通過(guò)頻譜兩點(diǎn)插值和迭代計(jì)算，實(shí)現(xiàn)了頻率高精度估計(jì)，后續(xù)稱為AM法。文獻(xiàn)[9]利用頻譜三點(diǎn)插值得到頻率估計(jì)值，后續(xù)稱為Candan法。針對(duì)單頻信號(hào)，AM法和Candan法分別是迭代類算法和非迭代類算法中綜合性能最好的算法，但處理多頻信號(hào)時(shí)，二者均受其他非待估計(jì)頻率分量頻譜泄漏的影響，頻率估計(jì)精度差，不能直接使用。

為抑制頻譜泄漏影響，文獻(xiàn)[10-11]在AM法的基礎(chǔ)上，采用頻譜相減策略實(shí)現(xiàn)了頻譜泄漏校正，提高了頻率估計(jì)精度，但在信號(hào)頻率間隔較近和中高信噪比條件下的頻率估計(jì)精度有待提高。

文獻(xiàn)[12]首先用改進(jìn)的AM法求取多頻信號(hào)各頻率分量的頻率粗估計(jì)值，然后構(gòu)造信號(hào)中非待估計(jì)的所有頻率分量，并通過(guò)相減策略以抑制頻譜泄漏的影響，與Ye等研究中的YA法具有相當(dāng)?shù)墓烙?jì)性能。但該算法減去非待估計(jì)頻率分量這一步驟只進(jìn)行了一次，導(dǎo)致頻率估計(jì)精度受頻率粗估計(jì)值的影響。特別是在信號(hào)頻率間隔較近、頻譜粗估計(jì)存在較大偏差時(shí)，由于信號(hào)中的非待估計(jì)頻率分量減去不徹底，使得算法仍受頻譜泄漏的影響，降低了頻率估計(jì)精度。

Djukanovi[13]通過(guò)濾除信號(hào)中非待估計(jì)頻率分量的方式來(lái)提高頻率估計(jì)精度，針對(duì)單頻實(shí)信號(hào)，提出了DFE-1法。首先利用Candan法對(duì)信號(hào)進(jìn)行預(yù)處理，得到頻率粗估計(jì)值，并生成參考信號(hào)；然后將信號(hào)與參考信號(hào)相乘以實(shí)現(xiàn)頻率搬移，通過(guò)濾除直流分量的形式降低頻譜泄漏的影響；最后采用AM法對(duì)濾除了非待估計(jì)頻率分量的信號(hào)進(jìn)行迭代計(jì)算，提高了信號(hào)頻率估計(jì)精度。針對(duì)多頻信號(hào)，在DFE-1法的基礎(chǔ)上進(jìn)行了推廣，提出了DFE-2法，實(shí)現(xiàn)了多頻信號(hào)頻率的高精度估計(jì)[14]。

DFE法對(duì)提高多頻信號(hào)的頻率估計(jì)精度提供了思路，但存在設(shè)計(jì)缺陷，影響了頻率估計(jì)精度。本文在分析頻譜泄漏影響和DFE法缺陷的基礎(chǔ)上，提出了新的頻率估計(jì)算法，并通過(guò)仿真實(shí)驗(yàn)驗(yàn)證了所提算法的有效性。

1 問(wèn)題描述

根據(jù)信號(hào)參數(shù)是否隨時(shí)間變化，采樣信號(hào)可分為平穩(wěn)信號(hào)和非平穩(wěn)信號(hào)，本文以平穩(wěn)信號(hào)為模型進(jìn)行分析。不失一般性，具有M個(gè)頻率分量的平穩(wěn)多頻信號(hào)如式(1)所示。

(1)

特別說(shuō)明：頻率分量個(gè)數(shù)M可以是已知的，也可以是未知的。在實(shí)際應(yīng)用中，M一般是未知的，可通過(guò)廣義阿卡克信息準(zhǔn)則和最小長(zhǎng)度描述法等方式計(jì)算得到。對(duì)頻率分量個(gè)數(shù)M進(jìn)行估計(jì)屬于信號(hào)檢測(cè)范疇，不屬于信號(hào)參數(shù)估計(jì)領(lǐng)域，是一個(gè)獨(dú)立的問(wèn)題。本文直接利用已有算法得到信號(hào)頻率分量個(gè)數(shù)，然后將M作為先驗(yàn)知識(shí)進(jìn)行處理，后續(xù)不再對(duì)其進(jìn)行討論。

特別地，當(dāng)M=2，且ω1=-ω2，a1=a2，θ1=θ2時(shí)，采樣信號(hào)稱為單頻實(shí)信號(hào)。

(2)

單頻實(shí)信號(hào)是一種特殊的多頻信號(hào)，含有正頻率分量和負(fù)頻率分量。設(shè)計(jì)多頻信號(hào)頻率估計(jì)算法時(shí)，對(duì)單頻實(shí)信號(hào)的頻率估計(jì)進(jìn)行討論，是有必要且具有代表性的。

在對(duì)多頻信號(hào)進(jìn)行頻譜分析時(shí)，可將其頻率表示為

(3)

式中：km=[ωmN/2π]為第m分量在頻譜中能量最大值點(diǎn)的索引； [·]為取最接近于·的整數(shù)；-0.5≤δm≤0.5為第m分量的頻譜偏移量。因此，要得到精確的頻率估計(jì)值，需要得到準(zhǔn)確的頻譜索引和精確的頻譜偏移量。

在對(duì)信號(hào)第m分量進(jìn)行頻率估計(jì)時(shí)，頻率估計(jì)精度受其他頻率分量頻譜泄漏和噪聲的疊加影響。由于DFT法具有很強(qiáng)的抗噪性，噪聲的影響相對(duì)較小。為直觀理解頻譜泄漏的影響，在無(wú)噪條件下，對(duì)多頻信號(hào)和單頻信號(hào)進(jìn)行頻譜分析，如圖1所示。

圖1 單頻信號(hào)和多頻信號(hào)頻譜

可以看出，受其他頻率分量頻譜泄漏的疊加影響，多頻信號(hào)待估計(jì)頻率分量的頻譜值大于同頻單頻信號(hào)的頻譜值，從而影響由頻譜分析得到的頻率估計(jì)值，與真實(shí)的頻率值存在估計(jì)偏差。

2 DFE算法

針對(duì)單頻實(shí)信號(hào)，Djukanovic提出了DFE-1算法，流程如表1所示。

表1 DFE-1算法流程

在此基礎(chǔ)上，Djukanovic提出了針對(duì)多頻信號(hào)的頻率估計(jì)算法，后續(xù)稱為DFE-2法。與DFT-1法相比較，兩個(gè)算法的流程一致。DFE-2法沿用了DFE-1法濾除信號(hào)中非待估計(jì)頻率分量的思路，但在進(jìn)行頻率粗估計(jì)和精估計(jì)時(shí)，均結(jié)合了Candan法和三點(diǎn)周期圖最大化法，沒(méi)再使用AM法。

DFE法抑制了非待估計(jì)頻率分量頻譜泄漏的影響，但存在設(shè)計(jì)缺陷：

(1) DFE-1法分別利用非迭代類Candan法和迭代類AM法對(duì)單頻實(shí)信號(hào)進(jìn)行頻率粗估計(jì)和精估計(jì)，但將粗估計(jì)和精估計(jì)兩個(gè)步驟完全分開(kāi)，既增加了計(jì)算量也降低了頻率估計(jì)精度。特別是在信號(hào)頻率較低，即信號(hào)正頻率和負(fù)頻率相隔很近、頻譜泄漏嚴(yán)重時(shí)，容易導(dǎo)致負(fù)頻率分量濾除不徹底，使得頻率估計(jì)精度受頻率粗估計(jì)的影響更加嚴(yán)重。

(2) 相位是關(guān)于頻率的函數(shù)，但DFE-1法在進(jìn)行頻率搬移時(shí)，只考慮了信號(hào)頻率，沒(méi)有考慮相位，降低了頻率估計(jì)精度。

(3) 針對(duì)多頻信號(hào)，DFE-2法分別進(jìn)行了一次頻率粗估計(jì)和頻率精估計(jì)，可以理解為進(jìn)行了迭代計(jì)算。但在處理濾除了非待估計(jì)頻率分量的單頻信號(hào)時(shí)，選用了非迭代類Candan法和三點(diǎn)周期圖最大化法，導(dǎo)致頻率估計(jì)精度受頻率粗估計(jì)的影響，且三點(diǎn)周期圖最大化法計(jì)算更加復(fù)雜，增加了計(jì)算量，降低了算法實(shí)時(shí)性。

(4) DFE-2法也沒(méi)有考慮相位的影響。

3 所提算法

為提高多頻信號(hào)的頻率估計(jì)精度，利用DFE法的思路，提出了新的參數(shù)估計(jì)算法，詳細(xì)步驟如下：

步驟1對(duì)采樣信號(hào)進(jìn)行FFT計(jì)算，并求信號(hào)頻譜中最大的M個(gè)極大值。

(4)

(5)

式中：f(·)為求函數(shù)·中最大M個(gè)極大值的索引；ki為第i分量的索引，i根據(jù)各頻率分量頻譜最大值遞減的順序依次排列，即k1和kM分別為信號(hào)頻譜M個(gè)極大值中最大值和最小值的索引。

(6)

(7)

式中，|·|和∠·分別為取復(fù)數(shù)·的模和角度。

步驟3對(duì)信號(hào)進(jìn)行頻率搬移，濾除信號(hào)中的其他非待估計(jì)頻率分量。

分析時(shí)，按照信號(hào)頻譜能量最大到最小的順序進(jìn)行分析，即從k1依次分析到kM。

首先構(gòu)造其他非待估計(jì)頻率分量的參考信號(hào)

(8)

式中，m根據(jù)非待估計(jì)分量頻譜最大值遞減的順序依次排列。

其次將參考信號(hào)與多頻信號(hào)相乘，將非待估計(jì)中的第m分量(m=1)搬移到0頻附近，得到搬移信號(hào)。

ym(n)=x(n)rm(n)

(9)

然后將搬移信號(hào)中0頻附近的信號(hào)能量視為直流分量，利用式(10)濾除直流分量，并將信號(hào)搬移回原頻率處，得到抑制了第m分量的降頻信號(hào)。

(10)

最后將降頻信號(hào)xM-1(n)代入式(9)和式(10)，濾除非待估計(jì)頻率分量中的第m分量(m=2)，得到降頻信號(hào)xM-2(n)。按照m的取值順序，重復(fù)式(9)、式(10)，依次濾除信號(hào)中非待估計(jì)頻率分量，最終得到只含有第i分量的降頻信號(hào)x1(n)。

步驟4采用AM法對(duì)降頻信號(hào)x1(n)進(jìn)行分析，在索引ki兩邊插值，間隔為0.5。利用式(11)計(jì)算插值點(diǎn)的頻譜值，并根據(jù)兩個(gè)插值點(diǎn)的頻譜由式(12)計(jì)算頻譜偏移量。

(11)

(12)

步驟5按照i的取值順序，循環(huán)計(jì)算步驟2～步驟4，得到每個(gè)頻率分量的參數(shù)粗估計(jì)值。

步驟6迭代計(jì)算步驟2～步驟5，進(jìn)一步提高各頻率分量的參數(shù)估計(jì)精度，得到各頻率分量的幅值和初相位估計(jì)值，并利用式(3)得到各頻率分量的頻率估計(jì)值。

綜上分析，算法的具體流程如表2所示。

表2 算法流程

此外，針對(duì)單頻實(shí)信號(hào)，根據(jù)表2的算法流程，即可得到濾除了負(fù)頻率分量頻譜泄漏影響的頻率估計(jì)值。

與DFE法的不同處在于：

(1) 不單獨(dú)區(qū)分信號(hào)頻率粗估計(jì)和精估計(jì)，均在迭代中進(jìn)行計(jì)算，既可以更加有效地濾除非待估計(jì)頻率分量的影響，提高頻率估計(jì)精度，也可以降低計(jì)算量，提升算法實(shí)時(shí)性。

(2) 濾除非待估計(jì)頻率分量時(shí)，考慮了相位的影響，有利于更加徹底地濾除非待估計(jì)頻率分量，從而提高頻率估計(jì)精度。

所提算法和DFE法、以及其他現(xiàn)有優(yōu)秀算法的性能將在仿真驗(yàn)證中進(jìn)行對(duì)比分析，檢驗(yàn)所提算法的優(yōu)越性。

4 仿真驗(yàn)證

為檢驗(yàn)所提算法的有效性，利用MATLAB軟件在不同條件下，對(duì)多頻信號(hào)和單頻實(shí)信號(hào)進(jìn)行頻率估計(jì)，且主要對(duì)多頻信號(hào)進(jìn)行分析。為降低計(jì)算時(shí)的隨機(jī)誤差，每組仿真進(jìn)行2 000次蒙特卡羅實(shí)驗(yàn)。為方便分析，將估計(jì)結(jié)果轉(zhuǎn)換為均方誤差(mean square errors, MSEs)，并用對(duì)數(shù)表示。

(13)

式中，L式蒙特卡羅實(shí)驗(yàn)次數(shù)。

實(shí)驗(yàn)時(shí)，以含有3個(gè)頻率分量的信號(hào)為例進(jìn)行頻率估計(jì)，并與YA法、AK法、DFE-1法、DFE-2法、Ye法[15]以及克拉美羅下限(Cramer-Rao lower bound，CRLB)[16]進(jìn)行對(duì)比分析。設(shè)采樣信號(hào)x(n)=1.5ej(3.1ωsn+θ1)+1.4ej((4.7+Δk)ωsn+θ2)+1.2ej((8.3+Δk)ωsn+θ3)+z(n)，長(zhǎng)度為128，索引間隔Δk以1為步長(zhǎng)從1增加到32，初相位θ1，θ2和θ3獨(dú)立隨機(jī)取值。

(1) 不同迭代次數(shù)

所提算法屬于迭代類算法，首先分析算法在不同迭代次數(shù)下的頻率估計(jì)性能。仿真時(shí)，設(shè)SNR為30 dB，迭代次數(shù)為1，2，3和4，結(jié)果如圖2所示。

圖2 不同迭代次數(shù)的頻率估計(jì)結(jié)果

經(jīng)由1次或2次迭代時(shí)，算法整體估計(jì)效果較差，但隨著頻率間隔增大而逐漸變好，且2次迭代的頻率估計(jì)精度高于1次迭代的頻率估計(jì)精度。當(dāng)頻率間隔較大(Δk≥5)時(shí)，3次和4次迭代具有相當(dāng)?shù)墓烙?jì)性能，信號(hào)各頻率分量頻譜間相互泄漏對(duì)所提算法的頻率估計(jì)精度影響非常小，頻率估計(jì)結(jié)果靠近CRLB。當(dāng)頻率間隔相距非常近時(shí)，4次迭代效果更好，考慮到算法的估計(jì)精度，特別是在信號(hào)頻率相隔較近、頻譜泄漏嚴(yán)重時(shí)頻率估計(jì)精度，后續(xù)實(shí)驗(yàn)均采用4次迭代計(jì)算。

(2) 無(wú)噪聲

由前文分析可知，利用頻域法對(duì)多頻信號(hào)進(jìn)行頻率估計(jì)時(shí)，受頻譜泄漏和噪聲的影響。因此，在無(wú)噪聲環(huán)境下，可以檢驗(yàn)各算法對(duì)頻譜泄漏的抑制能力，仿真結(jié)果如圖3所示。

圖3 無(wú)噪聲條件下的頻率估計(jì)結(jié)果

在不同頻譜泄漏程度下，AK法的頻率估計(jì)精度優(yōu)于YA法。當(dāng)頻率相隔較近、頻譜泄漏嚴(yán)重時(shí)，AK法優(yōu)于DFE-2法；當(dāng)頻率間的間隔較大、頻譜泄漏影響減小時(shí)，DFE-2法的效果優(yōu)于AK法。相比其他算法，本文算法采用頻率搬移策略濾除了信號(hào)中其他非待估計(jì)頻率分量，抑制了多頻信號(hào)中其他頻率頻譜間的泄漏影響，具有更強(qiáng)的頻譜泄漏抑制能力，明顯強(qiáng)于其他幾種算法，達(dá)到了算法的設(shè)計(jì)目的。

(3) 不同信噪比a

為說(shuō)明算法在不同信噪比下的估計(jì)性能，設(shè)x(n)=1.5ej(3.1ωsn+θ1)+1.35ej(7.2ωsn+θ2)+1.2ej(10.3ωsn+θ3)+z(n)，在SNR以1 dB為步長(zhǎng)從0增加到30 dB的條件下進(jìn)行了仿真實(shí)驗(yàn)，結(jié)果如圖4所示。

圖4 不同信噪比條件下的多頻信號(hào)頻率估計(jì)結(jié)果

當(dāng)信號(hào)頻率間隔小且各頻率分量的能量大小相近時(shí)，算法受噪聲影響大，在SNR<5 dB時(shí)的頻率估計(jì)精度很差。在5 dB15 dB、即頻率估計(jì)精度受頻譜泄漏影響更嚴(yán)重時(shí)，各算法均出現(xiàn)了頻率估計(jì)偏差，逐漸偏離CRLB，但相比于YA法、AK法和DFE-2法，本文算法的頻譜泄漏抑制能力更強(qiáng)，頻率估計(jì)精度最高，最接近CRLB，優(yōu)于其他幾種算法。

(4) 不同信噪比b

圖5 不同信噪比條件下的單頻信號(hào)頻率估計(jì)結(jié)果

當(dāng)信號(hào)頻率較低，即信號(hào)正頻率和負(fù)頻率分量相隔較近、頻譜泄漏嚴(yán)重時(shí)，本文算法優(yōu)勢(shì)明顯。當(dāng)SNR<5 dB時(shí)，Ye法、DFE-1法和本文算法具有相當(dāng)?shù)墓烙?jì)精度，均靠近CRLB。當(dāng)SNR>5 dB時(shí)，Ye法和DFE-1法的估計(jì)精度逐漸降低、逐漸偏離CRLB，且DFE-1法優(yōu)于Ye法。當(dāng)-5 dB

5 結(jié) 論

本文算法通過(guò)信號(hào)預(yù)處理、構(gòu)造參考信號(hào)、頻率搬移、濾除非待估計(jì)頻率分量等方式抑制了多頻信號(hào)中非待估計(jì)頻率分量頻譜泄漏的影響，并經(jīng)迭代計(jì)算得到了各頻率分量精確的頻率、幅值和初相位估計(jì)值。

仿真實(shí)驗(yàn)結(jié)果表明，所提算法有效地抑制了多頻信號(hào)中頻譜泄漏的影響，具有更高的頻率估計(jì)精度，頻率估計(jì)值的均方誤差比DFE法和其他現(xiàn)有優(yōu)秀算法的頻率估計(jì)結(jié)果更靠近CRLB。