基于DFT的相位差估計(jì)精度與改進(jìn)方法*

楊輝躍,涂亞慶,張海濤

(后勤工程學(xué)院信息工程系,重慶 401311)

?

基于DFT的相位差估計(jì)精度與改進(jìn)方法*

楊輝躍,涂亞慶*,張海濤

(后勤工程學(xué)院信息工程系,重慶 401311)

相位差是傳感器信號(hào)處理中重要的檢測(cè)參數(shù)。針對(duì)相位差高精度估計(jì)要求,在闡述DFT相位差估計(jì)原理基礎(chǔ)上,分析了影響估計(jì)精度的主要因素,推導(dǎo)出估計(jì)方差與信噪比、采樣長(zhǎng)度、頻率偏差及對(duì)稱窗型窗長(zhǎng)的具體關(guān)系,并給出了滿足精度要求的信噪比、采樣長(zhǎng)度和和頻率偏差條件。提出一種校正譜泄漏的相位差估計(jì)方法,先通過比值法計(jì)算出頻率偏差,然后考慮負(fù)頻率泄漏影響進(jìn)行相位差估計(jì),校正了短程和長(zhǎng)程兩類譜泄漏影響,給出了加矩形窗或Hanning窗的估計(jì)式和方法步驟。實(shí)驗(yàn)結(jié)果驗(yàn)證了估計(jì)精度分析及本文方法性能,科氏流量計(jì)應(yīng)用實(shí)驗(yàn)表明了方法的工程可行性和實(shí)用價(jià)值。

離散傅里葉變換;相位差估計(jì);頻譜泄漏;負(fù)頻率

相位差估計(jì)技術(shù)在傳感器信號(hào)處理、儀器儀表、故障診斷、電力電子等諸多領(lǐng)域有著廣泛且重要的應(yīng)用[1-5]。實(shí)際應(yīng)用中對(duì)相位差估計(jì)精度有較高的要求,例如,對(duì)稱負(fù)載下同步發(fā)電機(jī)三相輸出電壓之間的相位差在±0.6°范圍內(nèi),正常情況下電容型高壓電氣設(shè)備的介質(zhì)損耗角小于0.1°。因此,分析影響相位差估計(jì)精度的主要因素,實(shí)現(xiàn)對(duì)相位差高精度估計(jì)具有重要的現(xiàn)實(shí)意義和應(yīng)用價(jià)值。

目前針對(duì)相位差估計(jì),已提出多種估計(jì)方法,它們有各自的優(yōu)缺點(diǎn)和適用范圍。基于硬件電路的過零檢測(cè)法[6]根據(jù)兩路信號(hào)過零點(diǎn)時(shí)間計(jì)算相位差,計(jì)算速度快,但硬件成本高、抗干擾能力弱;數(shù)字相關(guān)法[7]利用兩路信號(hào)的相關(guān)函數(shù)估計(jì)相位差,對(duì)隨機(jī)噪聲抑制能力強(qiáng),但受諧波影響較大且要求整周期采樣;高階譜和互高階譜估計(jì)法[8]基于信號(hào)高階譜估計(jì)相位差,算法復(fù)雜,計(jì)算量大;Hilbert變換法[9]先對(duì)信號(hào)進(jìn)行90°移相,然后利用兩路信號(hào)移相前后的函數(shù)關(guān)系計(jì)算相位差,可動(dòng)態(tài)估計(jì)時(shí)變相位差,但易受諧波干擾。DFT相位差估計(jì)通過信號(hào)離散頻譜最大譜線處的相位相減獲得相位差,估計(jì)精度相對(duì)較高,可利用FFT快速算法,實(shí)時(shí)性強(qiáng)且便于硬件實(shí)現(xiàn),應(yīng)用廣泛。然而,DFT計(jì)算過程中的頻譜泄漏會(huì)影響相位差估計(jì)精度[10-12]。

本文在分析DFT相位差估計(jì)精度基礎(chǔ)上,提出一種校正譜泄漏的相位差估計(jì)方法。該方法首先通過比值法計(jì)算出頻率偏差,然后考慮負(fù)頻率泄漏影響,校正短程和長(zhǎng)程兩類譜泄漏[13-14]估計(jì)相位差,具有較高精度,其有效性和估計(jì)性能將在對(duì)比實(shí)驗(yàn)中進(jìn)行驗(yàn)證。

1 加窗DFT相位差估計(jì)原理

設(shè)2路單頻信號(hào)為s1(n)、s2(n),f0表示頻率,幅值和初始相位分別為A1、A2和θ1、θ2,以頻率fs(fs>2f0)對(duì)s1(n)和s2(n)進(jìn)行離散采樣,并用長(zhǎng)度為N的離散對(duì)稱窗w(n)進(jìn)行截短,得到有限長(zhǎng)離散加窗采樣序列:

s1w(n)=A1cos(2πf0n/fs+θ1)w(n)s2w(n)=A2cos(2πf0n/fs+θ2)w(n)

(1)

先估計(jì)相位θ1。對(duì)s1w(n)進(jìn)行DFT,得離散頻譜

(2)

忽略負(fù)頻率成分,在非整周期截?cái)鄷r(shí),即N≠mfs/f0,m∈Z+,頻譜泄漏存在頻率估計(jì)偏差δ,|δ|≤0.5。從而,可設(shè)f0=(k+δ)Δf,k為正整數(shù),Δf為頻率分辨率,則:

(3)

從而,s1w(n)的相位為:φ1=θ1+πδ。同理可得s2w(n)的相位:φ2=θ2+πδ。則估計(jì)相位差為:

(4)

2 估計(jì)精度分析

2.1 估計(jì)方差

在加性噪聲背景下,用長(zhǎng)度為N對(duì)稱窗wN(n)對(duì)含噪觀測(cè)信號(hào)x(n)進(jìn)行加窗截短,得:

xw(n)=x(n)wN(n)+z(n)wN(n)

(5)

對(duì)xw(k)進(jìn)行DFT,忽略負(fù)頻率成分,只考慮DFT頻譜的前N/2點(diǎn),有

(6)

(7)

利用二階泰勒公式展開,忽略高階無窮小,得:

(8)

設(shè)Pw為窗的平均功率,對(duì)噪聲譜有

(9)

由于高斯白噪聲任意兩條譜線的實(shí)部和虛部相互獨(dú)立,則

(10)

從而,根據(jù)信噪比定義和對(duì)稱窗頻譜,可得相位譜均值和方差

(11)

(12)

可見,相位差估計(jì)精度與信噪比、頻率偏差、所加對(duì)稱窗形狀和長(zhǎng)度有關(guān)。

2.2 滿足精度要求的估計(jì)條件

(13)

(14)

整理得:

(15)

上式即為滿足精度要求的DFT相位差估計(jì)條件。

(16)

由式(16)即可求得頻偏δ的取值范圍。

圖1 單頻信號(hào)負(fù)頻率泄漏

2.3 負(fù)頻率泄漏影響

上述相位差估計(jì)原理忽略了負(fù)頻率成分的影響,然而研究表明當(dāng)信號(hào)頻率很低或接近Nyquist頻率時(shí),由于頻譜泄漏導(dǎo)致負(fù)頻率譜峰疊加到正頻率譜,造成負(fù)頻率干涉影響相位差估計(jì)精度。如圖1所示,信號(hào)頻率在(0,π)的正頻率區(qū)間,(-π,0)和(π,2π)是正頻率區(qū)間關(guān)于ω=0和ω=π的對(duì)稱頻譜鏡像,稱為負(fù)頻率區(qū)間。當(dāng)信號(hào)頻率較低或接近Nyquist頻率時(shí),譜峰A靠近(0,π)兩端,旁瓣干涉迅速增大,甚至造成主瓣干涉,如圖1(b)和圖1(c)所示。此時(shí)考慮負(fù)頻率泄漏的影響十分必要。

3 校正方法

3.1 基本思想

根據(jù)上述分析,基于DFT的相位差估計(jì)精度與信噪比、頻率偏差、窗長(zhǎng)和窗形以及負(fù)頻率泄漏有關(guān)。因此,當(dāng)信噪比一定時(shí),可先進(jìn)行比值法獲得頻率偏差,然后不忽略負(fù)頻率成分計(jì)算相位差,從而實(shí)現(xiàn)相位差估計(jì)校正,改善估計(jì)精度。

3.2 方法原理

為獲得頻率偏差δ,本為采用比值法[15]即利用主瓣內(nèi)最大譜線和次大譜線的幅值譜之比進(jìn)行求解。設(shè)最大譜線和次大譜線分別為yk和yk+1,對(duì)應(yīng)相位為θk和θk+1,則有k

v=yk/yk+1=W(δ)/W(1+δ)

(17)

將所加窗函數(shù)頻譜代入式(14),利用峰值搜索或牛頓迭代法,即可求得頻率偏差δ。

若||θk-θk+1|-π|>ε,則認(rèn)為存在譜線干涉[16]。此時(shí),負(fù)頻率成分不可忽略。根據(jù)式(2),在不忽略負(fù)頻率的情況下重新推導(dǎo)相位差估計(jì)公式。若所加對(duì)稱窗為矩形窗,則經(jīng)推導(dǎo)后可得[17]

(18)

式中,φ1為S1(k0)的相位,c1=sin(2πk0/N);c2=sin[2π(k0+δ)/N];c3=2sin(πδ/N)sin[π(2k0+δ)/N]。

同理,對(duì)于第2路正弦采樣序列s2(n),有:

(19)

式中,φ2為S2(k0)的相位。由式(18)和(19),可求得相位差:

(20)

若所加窗為Hanning窗,同理可推導(dǎo)出校正后的相位差:

(21)

D1=1+cos(2π/N)-

2cos(πδ/N)cos[π(2k0+δ)/N]cos(2πk0/N)

D2=1+cos(2π/N)-

2cos(πδ/N)cos[π(2k0+δ)/N]cos[2π(k0+δ)/N]

3.3 方法步驟

方法實(shí)現(xiàn)步驟如下:

Step 1 加窗截短,用長(zhǎng)度為N對(duì)稱窗wN(n)對(duì)采樣信號(hào)x(n)進(jìn)行加窗截短;

Step 2 頻譜分析,對(duì)截取的信號(hào)進(jìn)行DFT計(jì)算,獲得信號(hào)離散頻譜;

Step 3 頻偏計(jì)算,根據(jù)信號(hào)離散頻譜,利用比值法,通過峰值搜索計(jì)算出頻偏δ;

Step 4 相位差估計(jì),根據(jù)窗wN(n)形狀選擇式對(duì)應(yīng)的(20)或式(21)估計(jì)相位差。

4 實(shí)驗(yàn)分析

采用單頻實(shí)正弦信號(hào)迭加高斯白噪聲,對(duì)本文理論推導(dǎo)進(jìn)行驗(yàn)證。兩路信號(hào)所加噪聲互不相關(guān),分別對(duì)頻偏、信噪比和采樣長(zhǎng)度、對(duì)稱窗型進(jìn)行單因素分析,并對(duì)比本文方法與原方法性能,給出精度要求下的估計(jì)條件。實(shí)驗(yàn)中,兩路單頻實(shí)正弦信號(hào)初始相位差Δθ=3.6°。

4.1 頻偏δ與窗函數(shù)的影響分析

在信號(hào)頻率f0=200 Hz,SNR=20 dB條件下,分別用矩形窗、Hanning窗、Hamming窗截取N=1 024點(diǎn)數(shù)據(jù),經(jīng)DFT估計(jì)信號(hào)相位差。同時(shí)采比值法校正后的相位差估計(jì)值作為比較,進(jìn)行200次獨(dú)立仿真實(shí)驗(yàn),得到相位差估計(jì)的均方根誤差與頻偏δ的關(guān)系如圖2所示。

圖2 頻偏δ與窗函數(shù)對(duì)估計(jì)精度的影響

圖2中3條線為按式(15)得出的均方誤差理論值,離散點(diǎn)為仿真結(jié)果。從圖2可以看出,仿真結(jié)果與理論值基本吻合。未經(jīng)校正的情況下,頻偏|δ|越接近零,即短程譜泄漏越小,DFT法相位差估計(jì)的均方根誤差越小。對(duì)比加不同窗的結(jié)果可知,加Hanning窗和加Hamming窗的均方誤差相差不大,|δ|較小時(shí),加矩形窗效果較好,當(dāng)|δ|較大時(shí),加Hanning窗和加Hamming窗的誤差更小。比值校正后,均方根誤差始終在某定值附近波動(dòng),與δ無關(guān)。該定值由信噪比和采樣長(zhǎng)度決定,在本實(shí)驗(yàn)給定條件下,該值約為0.253°。

4.2 信噪比與采樣長(zhǎng)度的影響分析

圖3為N=1 024,k0=200,δ=0.4條件下,校正后相位差估計(jì)均方根誤差與SNR的關(guān)系以及在k0=30,δ=0.4,SNR=20 dB條件下,校正后相位差估計(jì)均方根誤差與采樣序列長(zhǎng)度N的關(guān)系。從圖3可以看出,仿真結(jié)果與公式計(jì)算結(jié)果是吻合的,且SNR越高或N越大,相位差估計(jì)的均方根誤差越小。

圖3 信噪比SNR與采樣長(zhǎng)度對(duì)估計(jì)精度的影響

4.3 負(fù)頻率泄漏影響分析

考慮頻率較低時(shí)負(fù)頻率泄漏的影響。對(duì)2路單頻正弦信號(hào),分別用矩形窗和Hanning窗截取N=1 024點(diǎn),采樣頻率fs=1 000 Hz,則頻率分辨率Δf=fs/N=0.976 6 Hz。為更準(zhǔn)確地反映頻譜的內(nèi)在規(guī)律,以頻率分辨率為單位刻畫信號(hào)頻率,即f0在(0.5～2.5)Δf和(497.5～499.5)Δf內(nèi)取值,步長(zhǎng)為0.05Δf,仿真實(shí)驗(yàn)得到相位差估計(jì)的相對(duì)誤差與信號(hào)頻率的關(guān)系,如圖4所示。

圖4 相位差估計(jì)相對(duì)誤差

從圖中可以看出,式(20)和式(21)所示的負(fù)頻率修正DFT具有更高的估計(jì)精度,當(dāng)相對(duì)頻率f0/fd為整數(shù)時(shí),DFT法的誤差急劇下降,接近雙精度運(yùn)算的下限。這是因?yàn)楫?dāng)f0等于fd的整數(shù)倍時(shí),負(fù)頻率區(qū)間的譜峰A′的旁瓣在正頻率區(qū)間的譜峰A處的值正好為零,譜峰A并未受到任何影響。

4.4 滿足一定精度的估計(jì)條件計(jì)算

取α=0.01,查正態(tài)分布得u1-α/2=2.58,代入式(15)得誤差率以99%的概率落在20%內(nèi)的條件:

N·SNR≥166.41·Pw/|W(δ)|2

(22)

5 應(yīng)用驗(yàn)證

科氏流量計(jì)通過檢測(cè)兩路振動(dòng)信號(hào)的時(shí)間差(相位差)測(cè)量質(zhì)量流量,將本文方法用于科氏流量計(jì)進(jìn)行應(yīng)用驗(yàn)證。所用實(shí)驗(yàn)平臺(tái)如圖5所示,以10 kHz的采樣頻率,實(shí)際采集頻率約為146 Hz的流量計(jì)振動(dòng)信號(hào),分別用FFT法、比值校正DFT和本文方法(Hanning窗)進(jìn)行相位差估計(jì),每次計(jì)算64個(gè)采樣點(diǎn),所得結(jié)果換算成時(shí)間差如表1所示。其中理論參考值根據(jù)流量計(jì)性能曲線由實(shí)際流量換算而來。

圖5 科氏流量計(jì)實(shí)驗(yàn)平臺(tái)

表1 時(shí)間差估計(jì)值

對(duì)比可見,FFT法、比值校正DFT和本文方法估計(jì)精度依次提高,其中本文方法計(jì)算結(jié)果與理論參考值最為接近,驗(yàn)證了校正后的DFT相位差估計(jì)能有效抑制頻譜泄漏,較大幅度提高估計(jì)精度。

6 結(jié)束語

針對(duì)基于DFT的相位差估計(jì),推導(dǎo)出估計(jì)方差表達(dá)式,并論述了滿足一定精度要求的估計(jì)條件,在此基礎(chǔ)上,提出一種校正了頻譜泄漏的相位差估計(jì)方法,給出了方法原理和步驟,實(shí)驗(yàn)驗(yàn)證了估計(jì)精度分析的正確性以及本文方法的高精度性能,并在科氏流量計(jì)中進(jìn)行了應(yīng)用實(shí)驗(yàn)表明了本文方法可行性和價(jià)值。進(jìn)一步研究將基于本文方法,探討時(shí)變相位差的自適應(yīng)估計(jì)問題。

[1]Liu Xiao,Li Haisen,Zhou Tian.A 3-Subarray DOA Estimation Method Based on Phase Difference Measurement for Multi-Beam Bathymetry Sonar[J].Energy Procedia,2011(13):7683-7689.

[2]溫和,滕召勝,曾博,等.基于泄漏對(duì)消的電力諧波相角高精度估計(jì)算法[J].儀器儀表學(xué)報(bào),2009,30(11):2354-2360.

[3]楊輝躍,涂亞慶,張海濤,等.一種基于SVD和Hilbert變換的科氏流量計(jì)相位差測(cè)量方法[J].儀器儀表學(xué)報(bào),2012,33(9):2101-2107.

[4]伯恩,段發(fā)階,呂昌榮,等.基于CORDIC的交流相位跟蹤零差補(bǔ)償方法及其實(shí)現(xiàn)[J].傳感技術(shù)學(xué)報(bào),2014,27(2):198-203.

[5]馮文光,劉詩(shī)斌,李菊萍.數(shù)字磁通門傳感器的自動(dòng)相位對(duì)準(zhǔn)[J].傳感技術(shù)學(xué)報(bào),2012,25(2):212-214.

[6]龔國(guó)良,魯華祥.一種利用固定相移測(cè)量同頻正弦信號(hào)相位差的方法[J].儀器儀表學(xué)報(bào),2010,31(4):873-877.

[7]楊俊,武奇生,孫宏琪.基于相關(guān)法的相位差檢測(cè)方法在科氏質(zhì)量流量計(jì)中的應(yīng)用研究[J].傳感技術(shù)學(xué)報(bào),2007,20(1):138-145.

[8]寧輝.噪聲背景下正弦信號(hào)相位估計(jì)的互高階譜方法[D].長(zhǎng)春:吉林大學(xué),2001.

[10]張海濤,涂亞慶,牛鵬輝.相位差測(cè)量的FFT法和DTFT法誤差分析[J].電子測(cè)量與儀器學(xué)報(bào),2007,21(3):61-65.

[11]李炯,王巖飛.DFT相位估計(jì)算法及噪聲敏感頻率問題分析[J].電子與信息學(xué)報(bào),2009,31(9):2101-2103.

[12]Schuster S,Scheiblhofer S,Andreas S.The Influence of Windowing on Bias and Variance of DFT-Based Frequency and Phase Estimation[J].IEEE Transactions on Instrumentation and Measurement,2009,58(6):1975-1990.

[13]Young S S,Ho K J,Weui B W,et al.A Study on the Leakage Error in the Spectrum of Acoustic Intensity[J].JSME International Journal Series A,2004,47(1):42-46.

[14]刁瑞朋,孟慶豐.基于長(zhǎng)程泄漏補(bǔ)償?shù)牡逯礔FT方法及應(yīng)用[J].振動(dòng)與沖擊,2013,32(22):1-6.

[15]丁康,謝明,楊志堅(jiān).離散頻譜分析校正理論與技術(shù)[M].北京:科學(xué)出版社,2008.

[16]謝明,丁康,莫克斌.頻譜校正時(shí)譜線干涉的影響及判定方法[J].振動(dòng)工程學(xué)報(bào),1998,11(1):52-57.

[17]Tu Yaqing,Zhang Haitao,Mao Yuwen,et al.Unbiased Phase Delay Estimator with Negative Frequency Contribution for Real Sinusoids[J].Journal of Applied Sciences,2013,13(8):1160-1168.

Precision Analysis and Improvement Method for DFT Based Phase Difference Estimation*

YANGHuiyue,TUYaqing*,ZHANGHaitao

(Department of Information Engineer,Logistical Engineering University,Chongqing 401311,China)

Phase difference is an important parameter waiting for estimation in sensor signal processing.To fulfill precision requirements,effects against DFT based phase difference measurement method are analyzed on top of its measurement principle introduction.We derived the relationship between estimated variance with SNR,sampling length,frequency deviation,types and length of symmetrical window.Requirements for SNR,sampling length and frequency deviation are also put forward.In succession,a new phase difference measurement method with spectrum leakage revised is proposed.In the method,frequency deviation is firstly calculated by interpolation algorithm,and then phase difference is computed with the negative frequency influence considered.The short-range and long-range spectrum leakages are both taken into consideration in this method.Procedures and formulas for phase difference computation with rectangular window and Hanning window are given out.Experimental results prove the validity of precision analysis.Good performance in precision of the proposed method is also verified.Feasibility and practical value of the proposed method is shown by application experiments in Coriolis mass flow meter.

discrete Fourier transform;phase difference estimation;spectrum leakage;negative frequency

楊輝躍(1987-),男,湖南邵陽人,博士研究生,主要研究方向?yàn)橹悄軝z測(cè)與智能控制,數(shù)字信號(hào)處理,huiyue_yang@163.com;

涂亞慶(1963-),男,重慶人,教授,博士,博導(dǎo),主要研究方向?yàn)橹悄軝z測(cè)與智能控制,智能自動(dòng)化系統(tǒng),數(shù)字信號(hào)處理,yq.tu@163.com。

項(xiàng)目來源:國(guó)家自然科學(xué)基金項(xiàng)目(61271449,61302175);重慶市自然科學(xué)基金項(xiàng)目(CSTC2012jjA040006,CSTC2013jcyjA40030)

2014-10-08 修改日期:2014-11-15

C:7230

10.3969/j.issn.1004-1699.2015.01.017

TN911.7

A

1004-1699(2015)01-0093-06