CN105738696A - 全相位时移相位差频率估计方法及装置 - Google Patents

Abstract

本发明属于数字信号处理领域，为提供改进的全相位时移相位差频率估计方法及装置，能够提高精度、鲁棒性，降低信噪比门限。本发明采用的技术方案是，全相位时移相位差频率估计方法，包括如下步骤：Step 1输入M＝2N?1+L个样本即x(?N+1?L)～x(N?1)，获得第一路相位值Step 2对x(?N+1?L)～x(N?1?L)做全相位FFT，直接取其峰值谱位置索引值k＝k₂ ^*处的相位，作为第二路相位值Step 3将两路相位的差值乘以1/(△ωL)后，获得初始频偏估计Step 4对初始频偏估计进行频移补偿处理，得到修正后的频偏值估计Step 5输出归一化的频率估计值本发明主要应用于数字信号处理。

Description

全相位时移相位差频率估计方法及装置

技术领域

本发明属于数字信号处理领域，涉及给定有限样本，提升信号频率估计精度的问题。具体讲,涉及改进的全相位时移相位差频率估计法及其装置。

背景技术

频率估计是信号处理的最基本问题之一，也是振动分析[1]、阵列信号处理[2]、电磁学[3]等众多领域的关键问题。众所周知，该问题的最大似然估计落在信号的离散时间傅立叶变换(Discrete Time Fourier Transform,DTFT)的谱峰位置处，然而理想的DTFT工程上无法实现。故实际应用中，通常用DFT(Discrete Fourier Transform，其快速算法为FFT)取代DTFT。为克服FFT的栅栏效应，故需进一步对FFT峰值谱附近的谱线做内插[4,5]、校正[6]等处理，获得更精确的频率估计。

频率、幅值、相位是谐波信号的3个基本参数，而其中频率和相位的关系最为密切(频率对时间的积分即为相位变化量)。因而，从存在时延关系的两段样本中提取FFT峰值谱的相位差信息，就可以衍生出不同形式的相位差频率估计法[6-13]；而且基于相位差的频率估计已被理论证明是无偏的[13]，这也是相比于其他谱校正法的优势。故在有些对相位参数敏感的场合，如电力系统谐波分析[14]、光学工程[15]等，基于相位差的频率估计法往往起到其他频率估计法无法替代的作用。

全相位FFT(all-phase FFT,apFFT)因其兼备优良的抑制谱泄漏性能和很高的相位测量精度(即“相位不变性”[16])，故适合做基于相位差的频率估计。对存在时延关系的两段样本的apFFT峰值谱相位值作差分，即可得到全相位时移相位差频率估计法[9]，其价值已在密集谱校正[17]、光学目标重构[18]、电力谐波检测[14,19]等领域得到证实，文献[11]还推导出了全相位时移相位差法在时延为1个样本时的算法表述形式，提高了其实用性。为提升上述领域的检测性能，需对原估计器做进一步改进。

对于原有的全相位时移相位差频率估计器，文献[11]发现了其精度与信号的频率偏离程度有关(即频偏值δ越接近0.5，精度越低)，但该文献没有结合apFFT的性质对该实验结果做深入分析，也没有设计有效的改进措施。另外，衡量频率估计器精度的最有信服力的参考就是克拉美罗限(Cramer-Rao Lower Bound,CRB)[20]，文献[9,11]仅对全相位时移相位差法在噪声中误差做了测试，没有把频率估计方差与CRB进行对照。

针对以上不足，本发明提出一种高精度全相位时移相位差频率估计器，该估计器在全相位预处理环节引入无窗卷积窗替代文献[9,11]的双窗卷积窗；另外，对原有估计器的输出进行频移补偿，使得不管频率怎样偏离，峰值谱总能接近理想幅值，从而提高估计器抵御噪声的能力。改进后的估计器在任何频偏情况下，其精度不仅大大高于经典的全相位时移相位差估计器，而且高于国际上近期提出的Candan估计器[21]及其改进版本[22]，其频率估计方差在任意频偏情况下都可逼近CRB。

本发明提出的高精度全相位时移相位差频率估计法，在估计精度和噪声鲁棒性两方面大大改进了经典全相位时移相位差法的性能，在阵列信号处理、雷达通信、电磁学等领域有着非常广泛的应用前景。

参考文献

[1]Urbanek J.,Barszcz T.,Antoni J.A two-step procedure for estimationof instantaneous rotational speed with large fluctuations[J].MechanicalSystems and Signal Processing,2013,38(1):96-102.

[2]Wang W-Q,So H-C.Transmit subaperturing for range and angleestimation in frequency diverse array radar[J].Signal Processing,IEEETransactions on,2014,62(8):2000-11.

[3]Xu Y.,Yuan Q.,Zou J.,et al.Analysis of triangular periodic carrierfrequency modulation on reducing electromagnetic noise of permanent magnetsynchronous motor[J].Magnetics,IEEE Transactions on,2012,48(11):4424-7.

[4]齐国清,贾欣乐.插值FFT估计正弦信号频率的精度分析[J].电子学报,2004,32(4):625-9.Qi Guoqing,Jia Xinle.Accuracy Analysis of Frequency Estimation ofSinusoid Based on Interpolated FFT[J].Acta Electronica Sinica,2004,32(4):625-9.(In Chinese)

[5]Aboutanios E.,Mulgrew B..Iterative frequency estimation byinterpolation on Fourier coefficients[J].Signal Processing,IEEE Transactionson,2005,53(4):1237-42.

[6]丁康,钟舜聪.通用的离散频谱相位差校正方法[J].电子学报,2003,31(1):142-5.Ding Kang,Zhong Shuncong[J].A Universal Phase Diference CorrectingMethods on Discrete Spectrum,Acta Electronica Sinica,2003,31(1):142-5.(InChinese)

[7]齐国清,贾欣乐.基于DFT相位的正弦波频率和初相的高精度估计方法[J].电子学报,2001,29(9):1164-7.Qi Guoqing,Jia Xinle.High-Accuracy Frequency andPhase Estimation of Single-Tone Bas ed on Phase of DFT.Acta ElectronicaSinica,2001,29(9):1164-7.(In Chinese)

[8]齐国清.利用FFT相位差校正信号频率和初相估计的误差分析[J].数据采集与处理,2003,18(1):7-11.Qi Guoqing.Error Analysis of Frequency and PhaseEstimations Based on Phase Difference of Segmented FFTs[J].Journal of DataAcquisition&Processing,2003,18(1):7-11.(In Chinese)

[9]黄翔东,王兆华.全相位时移相位差频谱校正法[J].天津大学学报,2008,41(7):815-20.Huang Xiangdong,Wang Zhaohua.All-phase Time-shift Phase DifferenceCorrecting Spectrum Method[J].Transactions of Tianjin University,2008,41(7):815-20.(In Chinese)

[10]黄翔东,王兆华.基于全相位频谱分析的相位差频谱校正法[J].电子与信息学报,2008,30(2):293-7.Huang Xiangdong,Wang Zhaohua.Huang Xiangdong,WangZhaohua.All-phase Time-shift Phase Difference Correcting Spectrum Method[J].Journal of Electronics&Information Technology,2008,30(2):293-7.(In Chinese)

[11]张涛,任志良,陈光.改进的全相位时移相位差频谱分析算法[J].系统工程与电子技术,2011,33(7):1468-72.Zhang Tao,Ren Zhiliang,Chen Guang,et al.ImprovedSpectrum Analysis Algorithm for All-phaseTime-shift Phase Difference[J].Systems Engineering and Electronics,2011,33(7):1468-72.(In Chinese)

[12]安国臣,张秀清,王晓君.离散频谱的加权相位差校正方法[J].信号处理,2012,28(1):99-104.An Guochen,Zhang Xiuqing,Wang Xiaojun,et al.Weighting PhaseDifference Correcting Methods on Discrete Spectrum[J].Signal Processing,2012,28(1):99-104.(In Chinese)

[13]黄翔东,孟天伟,丁道贤.前后向子分段相位差频率估计法[J].物理学报,2014,63(21):214304-1-214304-7.Huang Xiangdong,Meng Tianwei,Ding Dao Xian,etal.A novel phase difference frequency estimator based on forward and backwardsub-segmenting[J].Acta Electronica Sinica,2014,63(21):214304-1-214304-7.(InChinese)

[14]杨宇祥,樊巨宝.全相位时移相位差法在电力谐波检测中的应用[J].电测与仪表,2012,49(7):24-28.Yang Yu Xiang,Fan Jubao.Application of the All-phaseTime-shift Phase Difference Method in the Power Harmonic Detection[J].Electrical Measurement&Instrumentation,2012,49(7):24-28.(In Chinese)

[15]杨颖,李醒飞,李洪宇.基于激光自混合效应的加速度传感器[J].光学学报,2013,33(2):234-40.Yang Ying,Li Xingfei,Li Hongyu,et al.Accelearation SensorBased on Laser Self-Mixing Interference[J].Acta Opitica Sinica,2013,33(2):234-40(In Chinese).

[16]王兆华,黄翔东.数字信号全相位谱分析与滤波技术[M].北京:电子工业出版社,2009.Wang Zhaohua,Huang Xiangdong.Techniques of All-phase Spectralanalysis and filtering fof digital signals[M].Beijing:Electronic IndustryPress,2009.(In Chinese)

[17]黄翔东,王兆华,罗蓬.全相位FFT密集谱识别与校正[J].电子学报,2011,39(1):172-7.Huang Xiangdong,Wang Zhaohua,et al.Discrimination and Correctionfor Dense All phase FFT Spectrums [J].Acta Electronica Sinica,2011,39(1):172-7.(In Chinese)

[18]于树海,王建立,董磊.基于全相位谱分析的傅里叶望远镜外场实验数据处理[J].光学精密工程,2012,20(10):2275-82.Yu Shuhai,Wang Jianli,Dong Lei,etal.Field experiment data processing of Fourier telescopy based on all phasspectrum analysis[J].Optics and Precision Engineering,2012,20(10):2275-82.

[19]彭祥华,周群,曹晓燕.一种高精度的电网谐波/间谐波检测的组合优化算法[J].电力系统保护与控制,2014,42(23):95-101.Peng Xiang-hua,Zhou Qun,Cao Xiao-yan.A high precision combinational optimization algorithm of power gridharmonic/inter-harmonic signal detection[J].Power System Protection andControl,2014,42(23):95-101.(In Chinese)

[20]Rife D,Boorstyn R.Single tone parameter estimation from discrete-time observations[J].Information Theory,IEEE Transactions on,1974,20(5):591-8.

[21]Candan C.A method for fine resolution frequency estimation fromthree DFT samples[J].Signal Processing Letters,IEEE,2011,18(6):351-4.

[22]Candan C.Analysis and Further Improvement of Fine ResolutionFrequency Estimation Method From Three DFT Samples[J].IEEE Signal ProcessingLetters,2013,20(9):913-6.

发明内容

为克服现有技术的不足，本发明旨在提供改进的全相位时移相位差频率估计方法及装置，能够提高精度、鲁棒性，降低信噪比门限。本发明采用的技术方案是，全相位时移相位差频率估计方法，包括如下步骤：

Step 1输入M＝2N-1+L个样本即x(-N+1-L)～x(N-1)，对其中的2N-1个样本即x(-N+1)～x(N-1)作无窗全相位FFT，取其峰值谱位置索引值k＝k^*处的相位将其减去时延相位k^*ΔωL＝k^*2πL/N后，再作取模及相位调整处理，获得第一路相位值其中，Δω＝2π/N，为角频率分辨率值；

Step 2对x(-N+1-L)～x(N-1-L)做全相位FFT，直接取其峰值谱位置索引值k＝k₂ ^*处的相位，作为第二路相位值

Step 3将两路相位的差值乘以1/(ΔωL)后，获得初始频偏估计

Step 4对初始频偏估计进行频移补偿处理，得到修正后的频偏值估计

Step 5输出归一化的频率估计值

上述步骤中，Step 1中的“取模及相位调整”的处理过程为：

(1)将相位差模除2π，得到相位值

(2)对做如下相位调整，而得到第1路最终相位输出值

上述步骤中，Step 4中的“频移补偿”的处理过程如下：

给定小阈值ξ，对初始频偏估计作如下分情况处理：

(1)若则构造序列-N+1-L≤n≤N-1对x(n)做频移，即

x(n)←x(n)s(n)，重复Step 1～Step 3的过程，得到新的频偏估计δ′，输出

(2)若不满足即直接输出

全相位时移相位差频率估计装置，包括：模数转换模块、时延模块、数字信号处理器DSP，待测频率信号x(t)进行模数转换处理，得到第一路输入信号x(n)；然后将第一路输入信号x(n)经过时延模块得到第二路输入信号x(n+L)，两路输入信号分别输送到数字信号处理器DSP，数字信号处理器DSP内部设置有全相位FFT预处理模块、相位提取模块、判别模块、频移补偿，数字信号处理器DSP主控器从I/O端口读取2路输入数据，进入内部RAM；然后依次由全相位FFT预处理模块、相位提取模块进行处理得到频偏估计，然后输出到判别模块判断频偏估计是否满足给定阈值，若不满足，判别模块将该频偏估计经频移补偿模块以及时延模块形成两路新的输入信号再输入到DSP进行处理，直到判别模块判断频偏估计满足给定阈值，从而得到最终结果。

全相位FFT预处理模块、相位提取模块用于执行下列步骤：

Step 1输入M＝2N-1+L个样本即x(-N+1-L)～x(N-1)，对其中的2N-1个样本即x(-N+1)～x(N-1)作无窗全相位FFT，取其峰值谱位置索引值k＝k^*处的相位将其减去时延相位k^*ΔωL＝k^*2πL/N后，再作取模及相位调整处理，获得第一路相位值

Step 2对x(-N+1-L)～x(N-1-L)做全相位FFT，直接取其峰值谱位置索引值k＝k₂ ^*处的相位，作为第二路相位值

Step 3将两路相位的差值乘以1/(ΔωL)后，获得初始频偏估计

Step 4对初始频偏估计进行频移补偿处理，得到修正后的频偏值估计

Step 5输出归一化的频率估计值

本发明的特点及有益效果是：

本发明提出的改进的全相位时移相位差频率估计方法，若应用于频率估计和实际工程领域，可产生如下有益效果：

第一、本发明改进的频率估计器可以提高原有的全相位时移相位差频率估计器的对噪声的鲁棒性，降低其信噪比门限。

对于原有的全相位时移相位差频率估计器，使用双窗全相位FFT对输入数据预处理，谱泄漏程度较严重。而本发明改进的全相位时移相位差频率估计器使用无窗全相位FFT代替原有的双窗全相位FFT，利用其频偏|δ|较小时的泄漏小、谱峰突出的特点，保证较高的频率估计精度。同时，在获得初始频偏估计后，本发明方法相比于原有的估计器，又对其进行频移补偿处理，进一步提高对噪声的鲁棒性。由于频移措施可以把信号绝大多数能量聚集在无窗apFFT的单根谱峰上，因而相比于双窗apFFT情况，抵御噪声能力大大增强，降低了信噪比门限值SNR_th。

第二、本发明改进的频率估计器还可以克服原有的全相位时移相位差频率估计器存在的整周相位模糊问题。

本发明改进的频率估计器在获得时延后第二路信号的谱峰相位后，对其取模处理并进行相位调整，有效保证了第二路观测相位在[-π,π]的范围内，克服了整周模糊问题，从而保证获得精确的频偏估计。

第三、本发明改进的频率估计器可以提高原有的全相位时移相位差频率估计器的精度。

本发明改进的频率估计器用无窗全相位FFT代替原有的双窗全相位FFT，并在原有估计器基础上引入频移补偿措施，该措施可以使得检测到的谱峰位置更接近真实值。通过把频率估计方差与CRB进行对照可以发现，本发明改进的估计器的RMSE曲线更接近CRB曲线，说明本发明估计器估计误差较小，估计精度较高。

附图说明：

图1全相位时移相位差频率估计流程。

图2全相位FFT谱分析的流程(N＝4)。

图3无窗apFFT和双窗apFFT在不同频偏值下的振幅谱。

图4三种方法的频率估计均方误差曲线。

图5本发明的硬件实施图。

图6 DSP内部程序流图。

具体实施方式

本发明提出的改进的全相位时移相位差频率估计器具有如下性能：

1)本发明频率估计器可以提高原有的全相位时移相位差频率估计器的精度，其频率估计方差更加接近克拉美罗限；

2)本发明频率估计器可以提高原有的全相位时移相位差频率估计器的对噪声的鲁棒性，降低其信噪比门限。

3)本发明频率估计器还可以克服原有的全相位时移相位差频率估计器存在的整周相位模糊问题。

1.估计器总测量流程

本发明提出的全相位时移相位差频率估计器的频率测量方法的流程如图1所示。

图1的处理步骤如下：

Step 1输入M＝2N-1+L个样本即x(-N+1-L)～x(N-1)，对其中的2N-1个样本即x(-N+1)～x(N-1)作无窗全相位FFT，取其峰值谱k＝k^*处的相位将其减去时延相位k^*ΔωL＝k^*2πL/N后，再作取模及相位调整处理，获得第一路相位值

Step 2对x(-N+1-L)～x(N-1-L)做全相位FFT，直接取其峰值谱位置索引值k＝k₂ ^*处的相位，作为第二路相位值

Step 3将两路相位的差值乘以1/(ΔωL)后，获得初始频偏估计

Step 4对初始频偏估计进行频移补偿处理，得到修正后的频偏值估计

Step 5输出得归一化的频率估计值

上述步骤中，Step 1中的“取模及相位调整”模块的处理过程为：

(1)将相位差模除2π，得到相位值

(2)对做如下相位调整，而得到第1路相位输出值

上述步骤中，Step 4中的“频移补偿”模块的处理过程如下。

给定小阈值ξ，对初始频偏估计作如下分情况处理：

(1)若则构造序列-N+1-L≤n≤N-1对x(n)做频移，即x(n)←x(n)s(n)，重复Step 1～Step 3的过程(即图1虚框中的处理过程)，得到新的频偏估计δ′，输出

(2)若不满足即直接输出

2.各处理步骤的详细过程

(1)用全相位FFT预处理数据并进行谱分析

图1中的apFFT谱分析流程图如图2所示。

图2中，首先对输入样本x(n)作全相位预处理：即用长为(2N-1)的卷积窗w_c对输入数据加权，然后将间隔为N的数据两两叠加(中间元素除外)而形成N个数据y(0),y(-1),…,y(-N+1)；其次，对y(0),y(-1),…,y(-N+1)做FFT得离散谱Y(k)，取其幅值即为振幅谱|Y(k)|，取其相角即为相位谱

图2中的w_c由长为N的“前窗”f和翻转的“后窗”b卷积得到的长为2N-1的卷积窗，即

w_c(n)＝f(n)*b(-n),n∈[-N+1,N-1] (1)

令R_N表示矩形窗，若f＝b＝R_N，则称图2为无窗apFFT；若f＝b≠R_N，则称双窗apFFT。前窗f的傅立叶变换表示为幅度和相位分离的一般形式为

F(jω)＝F_g(ω)e^j(N-1)ω/2 (2)

(2)峰值谱搜索得到的两路相位输出

对于输入序列ω₀＝2πf₀＝2π(k^*+δ)/N，k^*∈Z⁺，-0.5≤δ＜0.5，其传统FFT谱X(k)为

$X\left(k\right)=e^{j\[{\mathrm{\θ}}_0+\frac{N-1}{2}(k^{*}-k)\mathrm{\Δ}\mathrm{\ω}\]}{a}_0{F}_g\left(\right(k^{*}-k+\mathrm{\δ}\left)\mathrm{\Δ}\mathrm{\ω}\right),\mathrm{\Δ}\mathrm{\ω}=2\mathrm{\π}/N,k\∈\[0,N-1\]---\left(3\right)$

文献[16]指出，的apFFT幅度谱与传统FFT幅度谱的平方成正比，而相位谱与频率无关，等于输入中间样点的瞬时相位，则输入序列x(n)的apFFT谱Y(k)为

$Y\left(k\right)=e^{{\mathrm{j\θ}}_0}{a}_0{F_g}^2\left(\right(k^{*}-k+\mathrm{\δ}\left)\mathrm{\Δ}\mathrm{\ω}\right),k\∈\[0,N-1\]---\left(4\right)$

特别地，在峰值谱位置k＝k^*处，其apFFT谱为

$Y\left(k^{*}\right)=e^{{\mathrm{j\θ}}_0}{a}_0{F_g}^2\left(\mathrm{\δ}\mathrm{\Δ}\mathrm{\ω}\right)---\left(5\right)$

因而谱峰处的相位值为

即等于中心样点的瞬时相位，这种性质称为“相位不变性”。

对于延时L个样本后的序列由“相位不变性”，其谱峰处相位值等于中心样点x(-L)的瞬时相位，即

(3)取模及相位调整

取式(7)和式(6)的差值，即可得频率的理论估计式

实际需考虑观测相位与理论相位的整周模糊问题。不妨假定θ₀即为第1段apFFT的观测相位满足但时延后式(7)的理论相位θ₀-ω₀L很可能超出[-π,π]的范围，故可待定某整数即下式成立使得

对式(9)中可观测部分做2π取模处理得即

然后对做如下相位调整，

经式(11)调整后，联立式(9)，则可得存在某整数满足

把式(12)的合并成某整数可得

频偏估计可由下式求得

由于观测相位和调整后的相位均在[-π,π]的范围内，因而有

考虑到-0.5≤δ＜0.5，联立式(15)，可得(14)中的只能在{-1、0、1}中取值，则较正后的频偏估计为

$\widehat{\mathrm{\&delta;}}=\left\{\begin{array}{cc}-\widetilde{\mathrm{\&delta;}}+N/L,& \widetilde{\mathrm{\&delta;}}\&GreaterEqual;0.5\\ -\widetilde{\mathrm{\&delta;}},& -0.5\&le;\widetilde{\mathrm{\&delta;}}<0.5\\ -\widetilde{\mathrm{\&delta;}}+N/L,& \widetilde{\mathrm{\&delta;}}\&le;-0.5\end{array}\right.---\left(16\right)$

(4)频移补偿

给定小阈值ξ，对初始频偏估计作如下分情况处理

(1)若则构造序列-N+1-L≤n≤N-1，对x(n)做频移，即x(n)←x(n)s(n)，重复Step 1～Step 3的过程(即图1中虚框中的处理过程)，得到新的频偏估计δ′，输出

(2)若不满足则直接输出

注意上述两种分情况处理的结果，对不同符号的大频偏有自动调小的作用：即当检测为较大的正值时，则时域调制操作使得频谱朝负方向搬移，恰好搬到k^*Δf＝k^*/N附近；当检测为较大的负值时，则时域调制操作使得频谱朝正方向搬移，也恰好搬到k^*Δf附近。无论哪种情况，频移后的信号频偏都趋于0，以下实验将表明：只要信号频偏趋于0，则无窗apFFT谱峰总可以获得大的幅值。故其抵御噪声干扰的能力比频移前都会提升，保证了频率估计的高精度。

实验结果分析

实验1卷积窗和频偏对频率估计精度的影响

apFFT谱泄漏程度，既与频偏值δ有关，又与卷积窗类型有关。不妨假设输入复信号为(其中N＝32，ω₀＝(k^*+δ)Δω，k^*＝5)，幅值为1。图3给出了频偏δ取不同值情况下，无窗apFFT和双窗apFFT(其w_c由两个汉明窗卷积而得)的振幅谱，表1列出了谱峰幅值。

表1两种加窗情况下的apFFT谱峰幅值

从图3和表1可总结出两个规律：

(1)当频偏绝对值|δ|较小时(δ＝0,δ＝0.15)，无窗apFFT的旁谱泄漏较小，k＝5处的谱峰非常突出，其幅值接近于理想值1(δ＝0时，则不存在旁谱泄漏，谱峰幅度即为理想值1；δ＝0.15时，其幅值也高达0.928)。而双窗apFFT的谱峰左右各存在1根明显的旁谱，即使δ＝0时，主谱线幅值仍只有0.717。

(2)当频偏绝对值|δ|较大时(δ＝0.4)，无窗apFFT的旁谱泄漏变得严重，k＝5峰值谱线周围存在3根较大幅值的旁谱线；而双窗apFFT的k＝4处的旁谱幅值虽然有所减小，但k＝6处的旁谱幅值却继续增大，谱峰仍不突出，其幅值只有0.561。

因而，就双窗apFFT而言，无论频偏δ怎样取值，其谱峰幅值与理想值都存在较大差距，这必然会降低估计器抵御噪声的鲁棒性。相反，对于无窗apFFT，却可以尽可能利用|δ|较小时的泄漏小、谱峰突出的特点，保证较高的频率估计精度。

实验2频率估计器性能对比实验

设定输入信号为{x(n)＝2exp[j(5+δ)2π/N·n+π/3)]+w(n)，-N+1≤n≤N-1+L}，令N＝64，时延L＝7，w(n)为高斯复噪声。分别用经典全相位时移相位差法、改进的Candan内插法^[22]和本发明提出的高精度全相位时移相位差法(阈值ξ设定为0.1)做频率估计。设定其频偏值δ分别为0.1、0.2、0.3、0.4，对于每种频偏情况做1000次Monte-Carlo测频试验，并统计均方误差(root-mean-square error,RMSE)。图4(a)～图4(d)给出了各种RMSE曲线，并与克拉美罗限CRB做对照。需指出，文献^[20]给出了给定N个样本、信噪比为ρ(单位能量比值，SNR＝10lgρ)的频率估计的克拉美罗限为

$CRB=\frac{3}{2{\mathrm{\&pi;}}^2\mathrm{\&rho;}N(N^2-1)},---\left(17\right)$

由于本实验各种算法都耗费等长的M＝2N-1+L个样本，因而相应的CRB为

$CRB=\frac{3}{2{\mathrm{\&pi;}}^2\mathrm{\&rho;}(2N+L-2)(2N+L-1)(2N+L)}.---\left(18\right)$

从图4(a)～图4(d)，可总结如下规律：

(1)对于任意频偏情况，总体说来，本发明提出的高精度全相位时移相位差法(‘▽’标记)精度最高、改进的Candan估计器(‘□’标记)次之、经典全相位时移相位差法(‘○’标记)的精度最低。对于本发明方法，无论频率怎样偏移，RMSE曲线均几乎完全挨着CRB的开方曲线(实线标记)。这充分反映了本发明估计器通过选定无窗卷积窗和引入频移措施，在提升经典全相位时移估计器精度的效果是非常显著的，而且几乎已经做到了极限。

(2)本发明提出的改进的全相位时移相位差频率估计器，相比于经典全相位时移估计器，还降低了其低信噪比区域的信噪比门限值SNR_th(如表2所示)，特别是在频偏取值较大时(δ＝0.3和δ＝0.4)，其信噪比门限的改善尤为明显。这可从理论上得以解释：由于频移措施可以把信号绝大多数能量聚集在无窗apFFT的单根谱峰k^*Δf上，因而相比于双窗FFT情况，抵御噪声能力大大增强，降低了信噪比门限值SNR_th。

表2信噪比阈值SNR_th比较

(3)本发明方法在小频偏情况下，其信噪比门限与改进Candan估计器相当；但在大频偏情况下(δ＝0.3和δ＝0.4)，本发明方法的信噪比门限仍高于Candan方法，这可看作是为获得高精度性能所付出的代价。可从理论上解释：在给定M＝2N-1+L个样本情况下，Candan估计器做的是M点的FFT，同样功率的噪声平摊到2N-1+L根谱线上，因而每根谱线受到的噪声污染小。而apFFT在全相位预处理后，做的是N点的FFT，噪声功率只能平摊到N根谱线上，因而峰值谱线受到的噪声污染要比Candan估计器更大些。

下面就本发明涉及的硬件实现进行说明。

由图5，待测信号x(t)经过模数转换器A/D得到第一路输入信号x(n)，将x(n)经过时延模块得到第二路输入信号x(n+L)，将该两路信号同时输入数字信号处理器，在处理器内部经过全相位FFT预处理，频移补偿等算法处理，得到待测信号的频率估计值，最后由输出驱动及外设显示电路将频率估计结果实时显示。

其中，图4中的数字信号处理器为本发明硬件实施的核心器件，在信号频率估计过程中，完成以下主要功能：

(1)调用内部核心算法，对本发明提出的解析式进行构建，完成全相位FFT预处理，解决整周模糊和频移补偿等操作；

(2)利用得到的精确频偏值和输入的数据，完成待测高频信号的频率估计；

(3)反馈功能：根据实际应用要求对输入信号x(n)进行搬移，以得到精确的频偏估计，最终得到频率估计值；

(4)将本算法对待测信号频率的估计结果实时输出至驱动和显示模块。

需指出，由于采用了数字化的估计方法，因而决定了图5系统的复杂度、实时程度和稳定度的主要因素并不是图5中DSP器件的外围连接，而是DSP内部程序存储器所存储的核心频率估计算法。

DSP器件的内部程序流程如图6所示。

本发明将所提出的“改进的全相位时移相位差频率估计算法”这一核心算法植入DSP器件内，基于此算法完成高精度、高效的信号的频率估计。

图6流程分为如下几个步骤：

(1)对待测频率信号x(t)进行模数转换处理，得到第一路输入信号x(n)。然后将x(n)经过时延模块得到第二路输入信号x(n+L)，

(2)然后，DSP主控器从I/O端口读取2路输入数据，进入内部RAM；

(3)按图1中本发明的处理过程对2路输入数据进行全相位FFT预处理，相位提取，得到初始频偏估计；

(4)频移补偿。判断初始频偏估计是否满足给定阈值。若满足，直接将其作为最终频偏估计；若不满足，程序返回，进行频偏校正，得到最终频偏估计；直至测量结果符合工程要求。

(5)运行核心算法后，即可得到频率估计值；

(6)通过DSP的输出总线输出至外部的显示设备，将频率测量结果进行显示。

需指出，由于采用了DSP实现，使得整个频率估计操作变得更为灵活，可根据信号所包含的各种分量的具体情况，通过编程灵活改变算法的内部参数设置。

Claims (5)

1.一种全相位时移相位差频率估计方法，其特征是，包括如下步骤：

Step 1输入M＝2N-1+L个样本即x(-N+1-L)～x(N-1)，对其中的2N-1个样本即x(-N+1)～x(N-1)作无窗全相位FFT，取其峰值谱位置索引值k＝k^*处的相位将其减去时延相位k^*ΔωL＝k^*2πL/N后，再作取模及相位调整处理，获得第一路相位值其中，Δω＝2π/N，为角频率分辨率值；

Step 2对x(-N+1-L)～x(N-1-L)做全相位FFT，直接取其峰值谱位置索引值k＝k₂ ^*处的相位，作为第二路相位值

Step 3将两路相位的差值乘以1/(ΔωL)后，获得初始频偏估计

Step 4对初始频偏估计进行频移补偿处理，得到修正后的频偏值估计

Step 5输出归一化的频率估计值

2.如权利要求1所述的全相位时移相位差频率估计方法，其特征是，Step 1中的“取模及相位调整”的处理过程为：

(1)将相位差模除2π，得到相位值

(2)对做如下相位调整，而得到第1路最终相位输出值

3.如权利要求1所述的全相位时移相位差频率估计方法，其特征是，Step 4中的“频移补偿”的处理过程如下：

给定小阈值ξ，对初始频偏估计作如下分情况处理：

(1)若则构造序列-N+1-L≤n≤N-1对x(n)做频移，即x(n)←x(n)s(n)，重复Step 1～Step 3的过程，得到新的频偏估计δ'，输出

(2)若不满足即直接输出

4.一种全相位时移相位差频率估计装置，其特征是，包括：模数转换模块、时延模块、数字信号处理器DSP，待测频率信号x(t)进行模数转换处理，得到第一路输入信号x(n)；然后将第一路输入信号x(n)经过时延模块得到第二路输入信号x(n+L)，两路输入信号分别输送到数字信号处理器DSP，数字信号处理器DSP内部设置有全相位FFT预处理模块、相位提取模块、判别模块、频移补偿，数字信号处理器DSP主控器从I/O端口读取2路输入数据，进入内部RAM；然后依次由全相位FFT预处理模块、相位提取模块进行处理得到频偏估计，然后输出到判别模块判断频偏估计是否满足给定阈值，若不满足，判别模块将该频偏估计经频移补偿模块以及时延模块形成两路新的输入信号再输入到DSP进行处理，直到判别模块判断频偏估计满足给定阈值，从而得到最终结果。

5.如权利要求4所述的全相位时移相位差频率估计装置，其特征是，全相位FFT预处理模块、相位提取模块用于执行下列步骤：

Step 1输入M＝2N-1+L个样本即x(-N+1-L)～x(N-1)，对其中的2N-1个样本即x(-N+1)～x(N-1)作无窗全相位FFT，取其峰值谱位置索引值k＝k*处的相位将其减去时延相位k^*ΔωL＝k^*2πL/N后，再作取模及相位调整处理，获得第一路相位值其中，Δω＝2π/N，为角频率分辨率值；

Step 2对x(-N+1-L)～x(N-1-L)做全相位FFT，直接取其峰值谱位置索引值k＝k₂ ^*处的相位，作为第二路相位值

Step 3将两路相位的差值乘以1/(ΔωL)后，获得初始频偏估计

Step 4对初始频偏估计进行频移补偿处理，得到修正后的频偏值估计

Step 5输出归一化的频率估计值