正弦信号的直接FFT参数估计与相位差分法对比研究
李辉 王岩飞
①
（中国科学院电子学研究所 北京 100080） ②
（中国科学院研究生院 北京 100039）
摘 要：该文首先研究了直接基于FFT的正弦信号参数估计问题，揭示了频率与初相估计间的相互联系，并对相位差分法的估值误差公式进行了推导和仿真验证。二种算法的对比说明相位差分法运算量小，可以在不高的信噪比下获得彼此独立的高精度参数估值，因此更加有利于工程的实现。
关键词：直接FFT参数估计；对分迭代搜索；相位差分法；估值误差 中图分类号：TN911.72
①②
①
The Contrastive Study between Direct FFT and Phase Difference
in Parameter Estimation of Sinusoidal Signal
①②①
Li hui Wang Yan-fei
①
(Institute of Electronics, Chinese Academy of Sciences, Beijing 100080, China) ②
(Graduate School of the Chinese Academy of Sciences, Beijing 100039, China)
Abstract：This paper firstly studies parameter estimation problem directly based on FFT. Thereby, the relationship between the frequency and phase estimation is exposed. Subsequently, the estimation error formula of phase difference arithmetic are deduced and validated by computer simulation. The compare of two methods shows that phase difference arithmetic possesses lesser calculation quantity. Simultaneously, it can gain highly accurate, mutually independent parameter estimation under low SNR. So phase difference arithmetic is easy to realize in engineering field much more.
Key words: Direct FFT parameter estimation; Half-divided repetition search; Phase difference arithmetic; Estimation error
1 引言
正弦信号的参数（频率或初相）估计在雷达、声纳以及电子对抗等领域都有着极其广泛的应用。例如在雷达探测中，回波的频率和初相估计精度直接决定了被测物体的径向速度和距离测量精度。最直接的正弦信号频率估计就是在信号FFT之后，首先搜索到谱峰，再进一步运用对分法等迭代搜索算法以得到频率的精确估计。在此基础上将频率估值代入FFT的计
[1]
算式就可以得到对应的初相估计。这种方法的好处就是直接利用FFT的概念完成正弦信号参数的估值，无需进一步复杂的推导和证明，直观明了，但是在实际的工程应用中却不是最优的。这是因为在此种算法中频率估计的误差直接影响初相估计的精确性，而获得频率的精确估值就需要足够多的迭代次数，这样就往往不能满足系统对实时性的要求。相位差分法是在FFT粗测结果上的进一步校正，这种算法无需在频谱的最大和次大谱线间进行频率的搜索，只需对采样点分组后进行两次FFT就可以在不高的信噪比下获得精度相当高的频率和初
[2]
相估值，而且初相和频率的估计精度是彼此独立的，十分有利于工程的实现。本文首先分析了直接使用FFT进行参数估计时频率估计误差对初相估计的影响，并利用数值分析的方法
对实际工程中特定频点的参数估值问题进行了研究，定量地给出了频率和初相估值误差间的对应关系。文章第二部分对高斯白噪声环境中的相位差分参数估计算法进行了严密的数学推导，得出了估值误差与采样点数和信噪比之间严格的解析关系，利用MATLAB得到的仿真结果证明了理论推导的正确性。两种算法的对比说明相位差分法无需频率的迭代逼近，只需两次FFT就可以在不高的信噪比下获得彼此独立的高精度频率和初相估值，在运算量上具有很大的优势，尤其适用于实时性要求高的场合。
2 直接FFT参数估计研究
设以采样率fS?1?t采样后的点频信号可以表示为
n?0,1,2,??,N?1 s(n)?cos[2?f0n?t??]，（1）
?，对s(n)做FFT，搜索到谱峰后再进一步运用对分迭代等搜索算法得到频率的精确估值f0并代入FFT的计算式，得
N?1n?0?)?s(n)?exp[?j2?f?n?t]S(f?00
?n?t)?js(n)sin(2?f?n?t)??s(n)cos(2?f?00n?0n?0N?1N?1 （2）
所以初相?的估值
??arg?S(f?)??0??
?N?1?n?t)? s(n)cos(2?f（3） 0????1n?0?tan??N?1???n?t)?s(n)sin(2?f?0???n?0??f?l时，真实的初相 ?和f不相等，存在着一定的偏差，当f一般情况下f0000??arg?S(f0)?
?N?1?n?tl)?s(n)cos(2?f0??? （4） ?1n?0?tan??N?1???n?tl)?s(n)sin(2?f?0???n?0??。由于从??和?的表达式中很难直接得出频率估值误差??f，即l?1时，???可见当f00与初相估值误差之间明确的解析关系，所以只能用数值分析的方法对某一特定频点的参数估
计问题进行研究。
已知某在研的多频连续波测距系统的主侧音和第一匹配音分别为236MHz和220MHz，下面就以这两个频点为例，对直接利用FFT的参数估计问题使用数值分析的方法进行深入的研
?；究。基本的思路如下：首先对信号序列做FFT并搜索到离散谱线的最大值k0和次大值k0?为端点，运用对分迭代搜索算法，就可以进行频率的逼近搜索。算法的基本流以k0、k0[3]
程如图（1）所示，其中?f?fSN，为FFT的频率分辨率；??0，为设定的频率估值误差。
对信号序列做FFT?搜索到离散谱线的最大值 和次大值k0k0 令 ，ki?00ki??k0k?k?2i?i?1S(ki?f)?S(ki??f) 成立否？YesNoki?ki?ki??1?kiki?1? ，2ki?ki?ki??1?ki?ki?1? ，2No? 成立否？ki?1?f?ki?1?f??Yes??k?ff0i?1图（1）对分迭代搜索算法流程
在不同的信噪比下，利用MATLAB实现频率的迭代搜索和初相的计算
[4]
，当
、（2）所示。 ?f?fSN?500MHz1024时，迭代的估频误差和相应的估相误差如表（1）
估频误差(kHz) 信噪比SNR SNR=10dB SNR=6dB SNR=3dB 214.8 29.3 92.8 31.7 1.2 -80.2° 11.2° -33.8° -11.2° -1.1° -77.9° 9.3° -31.8° -11.0° 1.5° -77.6° 12.5° -36.5° -11.1° 1.8° 表（1）220MHz信号估频与估相误差
估频误差(kHz) 信噪比SNR SNR=10dB SNR=6dB SNR=3dB -160.2 83.9 -38.1 22.9 -7.5 7.7 58.0° -30.4° 15.0° -9.3° 2.3° -1.2° 57.7° -31.8° 13.0° -7.7° 2.1° 1.6° 62.2° -30.6° 16.1° -9.5° 1.8° -1.7° 表（2）236MHz信号估频与估相误差
从表（1）、（2）可以看出，当迭代的次数不够，频率偏差比较大时，信号的信噪比无论高低，初相估值都存在很大的偏差，不具备使用价值；经过几次迭代以后，频率估值逼近真实值，相应的初相估值也逼近真实值，此时误差主要由信噪比决定。
3 相位差分法研究
相位差分法的基本原理是将频率为f0，初相?的正弦信号
s(t)?exp[j(2?f0??)] （5） 以采样率fS进行采样，设信号的记录时间为T，采样点数为N，将得到的序列分为前后两个等长的子序列s1(n)和s2(n)，分别求出它们的FFT，得
[5]
S1(k)?Akexp(j?k)
（6） S2(k)?S1(k)exp(j?f0T)
其中Ak和?k分别为幅度项和相位项，有如下形式
[6]
Ak?sin[?(k?f0T/2)]sin[2?(k?f0T/2)/N]
（7） ?k???(1?2/N)(f0T/2?k)?
由式（7）可知幅度最大值处对应的离散频率为k0?[f0T/2]，用?1和?2表示S1(k)和S2(k)在最大谱线处的相位，由文献[7]可知，利用????2??1可以对粗测频率
??k2fs?k?f进行校准，偏差的估计为 fk00N?????f （8）f ?2?从而信号频率f0的估值
??f??f?f0k?
（9） ??????k0????f2???同时得到初相?的估值
3N?2N?2 ?1??2 （10）
2N2N在加性白噪声背景下，原N点的采样序列可以表示为
???r(n)?s(n)?z(n) （11）
r(n)的前N2点DFT可表示为
Rk?Akej?k?bej?z?Rkexp[j?R(k)]，k?0,1,2,?,N2?1 （12）
对于较大的DFT输出信噪比，由文献[6]可知在峰值谱线处的?R(k)可近似表示为
?R(k)k?k??k?00b sin(?z??k0) （13）
Ak0且bej?z为复高斯白噪声序列，有
E[bej?z]?0,var[（14） bej?z]?N?z2/2 所以可以推出
（15） var[bsin(?z??k)]?var[bcos(?z??k)]?N?z2/4
因此
N?z2b （16） var[?R(k0)]?var[?k0?sin(?z??k0)]?Ak04Ak2将式（7）代入式（15），有
?z21 （17） var[?R(k0)]??Nsinc2(??)2N?sinc2(??)?SNR?等于k0?f0T/2，表示实际信号频率和最大谱线对应频率之间的相对偏差，因此
??(?0.5,0.5)。
利用式（8），可得频率估计的方差
?)?var(f0var(?2??1)2 (?f)2(2?)