基于DFT的相位差估計精度與改進方法*

楊輝躍,涂亞慶,張海濤

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

?

基于DFT的相位差估計精度與改進方法*

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

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

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

離散傅里葉變換;相位差估計;頻譜泄漏;負頻率

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

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

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

1 加窗DFT相位差估計原理

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

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

(1)

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

(2)

忽略負頻率成分,在非整周期截斷時,即N≠mfs/f0,m∈Z+,頻譜泄漏存在頻率估計偏差δ,|δ|≤0.5。從而,可設f0=(k+δ)Δf,k為正整數,Δf為頻率分辨率,則:

(3)

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

(4)

2 估計精度分析

2.1 估計方差

在加性噪聲背景下,用長度為N對稱窗wN(n)對含噪觀測信號x(n)進行加窗截短,得:

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

(5)

對xw(k)進行DFT,忽略負頻率成分,只考慮DFT頻譜的前N/2點,有

(6)

(7)

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

(8)

設Pw為窗的平均功率,對噪聲譜有

(9)

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

(10)

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

(11)

(12)

可見,相位差估計精度與信噪比、頻率偏差、所加對稱窗形狀和長度有關。

2.2 滿足精度要求的估計條件

(13)

(14)

整理得:

(15)

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

(16)

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

圖1 單頻信號負頻率泄漏

2.3 負頻率泄漏影響

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

3 校正方法

3.1 基本思想

根據上述分析,基于DFT的相位差估計精度與信噪比、頻率偏差、窗長和窗形以及負頻率泄漏有關。因此,當信噪比一定時,可先進行比值法獲得頻率偏差,然后不忽略負頻率成分計算相位差,從而實現相位差估計校正,改善估計精度。

3.2 方法原理

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

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

(17)

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

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

(18)

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

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

(19)

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

(20)

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

(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 方法步驟

方法實現步驟如下:

Step 1 加窗截短,用長度為N對稱窗wN(n)對采樣信號x(n)進行加窗截短;

Step 2 頻譜分析,對截取的信號進行DFT計算,獲得信號離散頻譜;

Step 3 頻偏計算,根據信號離散頻譜,利用比值法,通過峰值搜索計算出頻偏δ;

Step 4 相位差估計,根據窗wN(n)形狀選擇式對應的(20)或式(21)估計相位差。

4 實驗分析

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

4.1 頻偏δ與窗函數的影響分析

在信號頻率f0=200 Hz,SNR=20 dB條件下,分別用矩形窗、Hanning窗、Hamming窗截取N=1 024點數據,經DFT估計信號相位差。同時采比值法校正后的相位差估計值作為比較,進行200次獨立仿真實驗,得到相位差估計的均方根誤差與頻偏δ的關系如圖2所示。

圖2 頻偏δ與窗函數對估計精度的影響

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

4.2 信噪比與采樣長度的影響分析

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

圖3 信噪比SNR與采樣長度對估計精度的影響

4.3 負頻率泄漏影響分析

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

圖4 相位差估計相對誤差

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

4.4 滿足一定精度的估計條件計算

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

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

(22)

5 應用驗證

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

圖5 科氏流量計實驗平臺

表1 時間差估計值

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

6 結束語

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

[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]溫和,滕召勝,曾博,等.基于泄漏對消的電力諧波相角高精度估計算法[J].儀器儀表學報,2009,30(11):2354-2360.

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

[4]伯恩,段發階,呂昌榮,等.基于CORDIC的交流相位跟蹤零差補償方法及其實現[J].傳感技術學報,2014,27(2):198-203.

[5]馮文光,劉詩斌,李菊萍.數字磁通門傳感器的自動相位對準[J].傳感技術學報,2012,25(2):212-214.

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

[7]楊俊,武奇生,孫宏琪.基于相關法的相位差檢測方法在科氏質量流量計中的應用研究[J].傳感技術學報,2007,20(1):138-145.

[8]寧輝.噪聲背景下正弦信號相位估計的互高階譜方法[D].長春:吉林大學,2001.

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

[11]李炯,王巖飛.DFT相位估計算法及噪聲敏感頻率問題分析[J].電子與信息學報,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]刁瑞朋,孟慶豐.基于長程泄漏補償的迭代插值FFT方法及應用[J].振動與沖擊,2013,32(22):1-6.

[15]丁康,謝明,楊志堅.離散頻譜分析校正理論與技術[M].北京:科學出版社,2008.

[16]謝明,丁康,莫克斌.頻譜校正時譜線干涉的影響及判定方法[J].振動工程學報,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-),男,湖南邵陽人,博士研究生,主要研究方向為智能檢測與智能控制,數字信號處理,huiyue_yang@163.com;

涂亞慶(1963-),男,重慶人,教授,博士,博導,主要研究方向為智能檢測與智能控制,智能自動化系統,數字信號處理,yq.tu@163.com。

項目來源:國家自然科學基金項目(61271449,61302175);重慶市自然科學基金項目(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