游行远 杨平 徐彬彬

（1.哈尔滨工程大学 信息与通信工程学院， 黑龙江 哈尔滨150001；2.武汉船舶通信研究所， 湖北 武汉430079）

在短波突发通信系统中，发射机与接收机两端电台振荡器不稳定［1］，导致发送端与接收端载波之间存在偏差.短波传输环境中信道衰落导致的多普勒频移，在经过相干解调下变频之后接收到的基带信号存在频率偏差.现有短波调解器的波形格式参照了北约STANAG4285、4529、5066 标准，一般通过在波形格式起始位置或数据块结构之间插入已知的8PSK 符号序列（训练序列）来进行信号同步捕获、信道参数估计、自适应均衡.同时波形负荷的数据部分会根据信息速率的需求采用不同的调制方式，如QPSK、8PSK、16QAM、64QAM 等.

频率偏移会直接导致接收信号中有用的信号成分功率衰减，特别对PSK、QAM 调制来说，频偏会导致数据符号相位旋转与扩散，影响数据信息正确解调，增大系统误码率［2］.因此，需要在接收端对频偏进行有效的估计，以修正其对接收信号的影响，提高接收机的性能.

目前国内外已经提出了不少频率估计方法，时域方法中基于自相关的方法一直广受关注，其利用信号的自相关相位来提取频率.其中，文献［3］中通过将不同延迟的自相关相位进行加权平均来提取频率，但算法的频率估计范围受限于参与计算的最大自相关延迟.文献［4］中利用自相关相位的差分来实现频率估计，避免了频率受限问题，但低信噪比时性能下降，计算复杂度增加.在文献［4］算法的基础之上，文献［5］中通过对自相关函数进行归一化来改进算法，减少了算法求自相关相位的次数.频域方法大多是基于离散傅里叶变换（DFT）的方法，此类方法物理意义明确，计算量小，因而得到了广泛的应用［6］.现有基于DFT 的频偏估计方法主要分为粗估计与精估计两个步骤，粗估计利用DFT 对信号频谱峰值进行定位，精估计是目前此类方法的研究热点.如文献［7-9］中利用DFT 结果峰值谱线与相邻谱线进行一次插值来得到估计结果.此类插值方法简单易行，但性能相对较差.文献［10-11］中提出了利用DFT 结果进行二分搜索的方法，但需要在粗估计时通过补零的方式（即增加DFT 点数）来获取更精细的频率分辨率［12］，导致计算复杂度大幅提升.文献［13］中利用DFT 峰值谱线±1/2 刻度的傅里叶系数进行迭代，算法理论均方误差（MSE）为Cramer-Rao下界（CRLB）的1.0147 倍.

为适应不同速率调解器的波形格式，保证不同长度训练序列的频偏估计精度，文中提出了一种基于傅里叶系数插值的迭代频偏估计算法，首先进行非补零DFT 的粗估计（即DFT 点数与训练序列长度保持一致），同时计算信号的信噪比，然后根据信噪比信息、训练序列长度、对插值算法的误差分析来确定频偏搜索区间，最后在搜索区间内进行迭代，并利用抛物线插值获取更精确的结果.

1 系统模型

假设长度为L 的8PSK 符号训练序列为ck（k=1，2，…，L），理想同步条件下，经过匹配滤波后，接收到的采样信号可以表示为

式中：A 为信号幅度；Ts为符号采样间隔；f0为频率偏移；φ 为相位偏移，一般由信道时延与采样偏差产生；n（k）为高斯白噪声（均值为0，方差为σ2）.接收信号与本地训练序列逐点共轭相乘，可得

式中，（·）*表示求共轭，w（k）=n（k）c*k.不失一般性，设ckc*k=1，同时将频偏进行归一化，得到F=f0Ts（以下不做特殊描述时均指归一化频偏），式（2）可简化输出信号为

式中，w（k）与n（k）有相同的统计特性，服从高斯分布.因此，在短波突发通信中的频偏估计可以等效于高斯白噪声条件下复正弦信号的频率估计［1］.

在粗估计步骤中，对式（3）进行N 点DFT 变换（N≥L），可得

式中，W（n）为噪声的傅里叶系数（服从均值为0、方差为Nσ2的高斯分布）［14］，DFT 的归一化分辨率为1/N.在实际系统中，真实频偏一般表示为

式中：m 为整数部分，即DFT 峰值谱线的位置；δ∈［-1/2，1/2］，为小数部分，即真实频偏相对于DFT 峰值谱线的偏差.文中的目标是低信噪比条件下对δ 的精确估计.

同时利用接收信号与本地训练序列之间的相关性，可得

由式（6）可求得接收信号的信噪比:

2 算法原理

2.1 搜索区间的确定

根据频偏粗估计DFT 结果，其峰值谱线与相邻谱线可表示为

式中，p=-1，0，1，分别对应谱线的位置.为了方便描述，定义Yp=Y（m+p），Wp=W（m+p）.

当N≫p 时，

将式（9）代入式（8），可得

式中，B=Aexp（jφ） [1-exp（j2 δ） ],通过DFT 峰值谱线与相邻谱线进行Jacobsen 插值［8］，得

式中，WN=W-1-W1，WD=2W0-W-1-W1.高斯白噪声的傅里叶系数Wp为，而BN 为O（N），其中O（·）表示函数阶符号［13］.因此，当N 取较大值时，式（11）可近似于

由于噪声均值E［Wp］=0，故为δ 的无偏估计.由于噪声傅里叶系数的方差var［Wp］=Nσ2，噪声的虚部方差为var［Im（Wp）］=Nσ2/2，可得插值结果的近似方差为

式中：Q（·）为高斯Q 函数；△＞0，为搜索区间两端到δ 的距离.

置信度可以反映置信区间的可靠程度，置信度越大，用以估计实际频偏值出错的概率越小.然后在实际的应用中，除了置信区间的可靠以外，还需要保证足够的精确.文中取置信度为1-a，a=10-5，得△=，则δ 的置信区间为，将其作为精频偏估计的搜索区间，得到初始状态为

2.2 精频偏估计

精频偏估计主要利用迭代方法在区间［m0-△0，m0+△0］内进行搜索.

首先利用搜索区间的中心点及两端，求得对应的Y（n）（n=m0-△0，m0，m0+△0），根据这3 点的大小关系进行如下判决:

因此，下一次迭代的初始状态为m1=m0+δ0△0，△1=△0/4.

（2）当两端S（m0-△0）或S（m0+△0）为最大值时，真实频偏属于区间［m0-△0，m0-△0/2］或［m0+△0/2，m0+△0］，即，为了方便算法迭代，文中取δ0=-3/4或δ0=3/4，则下一次迭代的初始状态为m1=m0±3△0/4，△1=△0/4.

利用上述判决方法进行I 次迭代，可得精频偏估计结果：

算法的具体步骤如下：

1）粗频偏估计Y，并确定频谱峰值m，即同时根据式（7）计算信噪比.

2）确定搜索区间.利用式（11）对频谱峰值谱线与相邻谱线进行插值，得到，并利用插值算法误差分析的结果（式（13））来得到插值结果的近似方差，最后利用△=4.42来确定精频偏估计的搜索区间

3）精频偏估计.利用I 次（i=0，1，...，I-1）迭代插值，根据式（4）计算并根据其大小关系进行判决.

本次迭代结果为mi+1=mi+δi△i，△i+1=△i/4.I 次迭代结束后，最终结果为

3 算法性能仿真分析

为进一步验证所提算法的性能，文中对其进行仿真分析，仿真中训练序列采用8PSK 调制，加性高斯白噪声（AWGN）信道条件，粗频偏估计采用非补零的DFT，即N=L，通过Nmc=106次蒙特卡洛实验得到MSE，计算公式如下：

因短波通信中的频偏估计可以等效为高斯白噪声条件下复正弦信号的频率估计，故可以利用其CRLB 对算法的性能进行验证，CRLB 可表示为［16］

对于不同长度的训练序列，文中算法的MSE 及CRLB 理论线如图1所示，其中设定的δ 均匀分布于区间［-1/2，1/2］.从图中可知，不同长度的训练序列存在各自的SNR 阈值：L为128、256、512 的训练序列的SNR 阈值分别为-4、-6、-9dB 左右，即粗估计时谱线峰值准确定位的概率下降导致了低于SNR 阈值时性能急剧恶化，而高于SNR 阈值时性能趋于平稳［17］.当I=0，即未使用迭代算法仅使用一次插值获取的结果时，算法的MSE 无法逼近CRLB，而使用I=1，2 的迭代插值后，算法对不同长度的训练序列均可逼近CRLB.

当L=512、δ 取不同值时文中算法的MSE 如图2所示.当δ 为0.05、0.20、0.40 时，SNR 阈值分别为-11、-10、-9dB；当δ 较小时，算法能获得更小的阈值.

图1 文中算法对不同长度的训练序列的MSEFig.1 MSE of the proposed algorithm for training sequence with different length

图2 L=512、δ 取不同值时文中算法的MSEFig.2 MSE of the proposed algorithm with different values of δ and L=512

当L=512 时文中算法与现有算法的MSE 比较见图3，其中设定的δ 均匀分布于区间［-1/2，1/2］.时域算法集中利用不同延迟的自相关相位arg［Ryy（m）］来提取频率，其中arg［·］表示取幅角，Ryy（m）表示为

式中，M 为最大延迟.文献［3］算法的估计范围存在局限性，其估计范围为0≤≤1/M，文中取M=20.文献［4］算法解决了估计范围的问题，其M 需小于或等于L/2，为了提高算法性能，文中取M 为128 和256作比较.文献［5］算法先归一化Ryy（m）得到，再利用与进行差分相乘，并对其结果进行加窗后，将自相关相位作为估计结果.仿真结果表明：在时域类方法中，文献［3］算法在信噪比ρ∈［-20，10］时皆可对频率进行跟踪，但MSE 性能较差；文献［4］算法的MSE 性能在ρ∈［-20，-10］时较差，但随着SNR 的增加性能逐步改善，当M 为128和256 时，分别在SNR 为10、-7 dB 处逼近CRLB.当M=256 时文献［5］算法在SNR 为-5 dB 处逼近CRLB，其性能趋近于文献［4］算法，但因仅需要一次求相位运算而降低了计算复杂度.频域类方法中，文献［3］算法无法逼近CRLB;文献［10-11］算法在粗估计使用N=2L 的DFT 时，SNR 阈值为-11 dB，表明适当地增加DFT 长度可以提高DFT 频谱的精细程度，提高频谱峰值准确定位的概率，但会提升计算复杂度；文献［13］算法在I=2 时可以逼近CRLB.

图3 L=512 时不同算法的MSE 比较Fig.3 Comparison of MSE among different algorithms with L=512

当L=512、SNR 为0dB 时文中算法与现有算法的MSE 比较见图4，其中MSE 与CRLB 之比R 为

SNR 为0dB 时CRLB 的理论值为CRLB0.时域类方法在δ 取不同值时性能较稳定，而频域类方法在δ 取不同值时性能会起伏.其中文献［10-11］算法在使用非补零DFT 时无法逼近CRLB，而在使用N=2L 的DFT 时，通过10 次迭代搜索，针对不同的δ 均可以获得较好的性能.文献［13］算法在I=1、δ 越趋近于0时性能越好，而在I=2 时可逼近CRLB;文中算法在I=1，2 时皆可逼近CRLB.

SNR 为0 dB 时文中算法与文献［13］算法在不同插值次数（I′）下的MSE 与CRLB 之比R 如表2所示.由于文献［13］算法未使用粗估计的DFT 结果进行插值，仅利用其对频偏进行定位，故以插值次数为对比标准，相应地文中算法的迭代次数I=I′-1，文献［13］算法的I=I′.从表中可知，文中算法在δ=0.0，0.1，0.2 时通过两次插值即可获得较稳定的性能，得到的R 可达1.005 左右，而δ 较大时需要3 次插值才能获得较稳定的性能.文献［13］算法在δ 较小时需要3 次插值才能获得趋于稳定的性能，而δ较大时需要4 次插值才能获得趋于稳定的性能，其R 在1.01～1.02 之间.

图4 L=512、SNR=0 dB、δ 取不同值时几种算法的MSE比较Fig.4 Comparison of MSE among several algorithms with L=512，SNR=0dB and different values of δ

计算复杂度也是评价算法性能的重要指标［18］，时域类方法主要集中于自相关性Ryy（m）的计算，其计算复杂度与最大延迟有关，需要M（2L-M-1）/2 次复数相乘运算，当M=L/2 时，其计算复杂度为3L2/8，为O（L2）.频域类方法包括粗频偏估计与精频偏估计两部分，粗频偏估计主要为N=L 的DFT，需要进行（L/2）log2L 次复数运算；而精频偏估计主要利用傅里叶系数的相互关系进行迭代，主要集中于Y（n）的计算，计算量为CY=L（8CM+6CP），为O（L），CM为复数乘法运算，CP为复数加法运算，插值部分的计算复杂度多为O（1），可以忽略.文中算法通常迭代两次即可获得较好的性能，单次迭代计算量为3CY.文献［13］算法需要迭代4 次才能保证算法的性能趋于稳定，单次迭代计算量为2CY.文献［10-11］算法在粗估计时的补零DFT 需要Llog2（2L）次复数运算，同时其精频偏估计需要迭代10 次，单次迭代计算量为CY.

表2 两种算法的R 值比较Table2 Comparison of R between two algorithms

4 结论

通过对现有基于DFT 类方法的研究，文中分析了Jacobsen 插值算法的误差，建立了估计值方差与信噪比及序列长度的一种函数关系，据此来提取估计值的置信区间并作为精频偏搜索区间，之后利用迭代在区间内进行插值的精频偏估计.仿真结果表明：在确立的置信区间内进行迭代精频偏估计，能有效地提高算法的估计精度，减少迭代次数；在低信噪比条件下，通过两次迭代文中算法的MSE 即可逼近CRLB 理论线， 故适用于短波突发通信的信道环境.文中算法也适用于高斯白噪声环境下对复正弦信号的频率估计.

［1］Palmer J，Rice M.Low-complexity frequency estimation usingmultipledisjointpilotblocksinburst-modecommunications［J］.IEEE Transactions on Communications，2011，59（11）：3135-3145.

［2］Pollett T，Moeneclaey M.The effect of carrier frequency offset on the performance of band limited single carrier and OFDM signals ［C］//Proceedings of IEEE Global Telecommunications Conference.London：the Key to Global Prosperity，1996：719-723.

［3］Fitz M P.Further results in the fast estimation of a single frequency ［J］.IEEE Transactions on Communications，1994，42（2/3/4）：862-864.

［4］Mengali U,Morelli M.Data-aided frequency estimation for burst digital transmission ［J］.IEEE Transactions on Communications,1997,45（1）：23-25.

［5］Awoseyila A B,Kasparis C,Evans B G.Improved single frequency estimation with wide acquisition range［J］.Electronics Letters，2008，44（3）：245-247.

［6］曹燕.含噪实信号频率估计算法研究［D］.广州：华南理工大学电子与信息学院，2012.

［7］Quinn B G.Estimating frequency by interpolation using Fourier coefficients ［J］.IEEE Transactions on Signal Processing，1994，42（5）：1264-1268.

［8］Jacobsen E，Kootsookos P.Fast，accurate frequency estimators ［J］.IEEE Signal Processing Magazine，2007，24（3）：123-125.

［9］CandanC.Amethodforfineresolutionfrequencyestimation from three DFT samples ［J］.IEEE Signal Processing Letters，2011，18（6）：351-354.

［10］Zakharov Y V，Tozer T C.Frequency estimator with dichotomous search of periodogram peak［J］.Electronics Letters，1999，35（19）：1608-1609.

［11］Zakharov Y V，Baronkin V M，Tozer T C.DFT-based frequency estimators with narrow acquisition range ［J］IEE Proceedings：Communications，2001，148（1）：1-7.

［12］Fang L，Duan D，Yang L.A new DFT-based frequency estimator for single-tone complex sinusoidal signals ［C］//Proceedings of Military Communications Conference.Orlando：IEEE，2012：1-6.

［13］Aboutanios E，Mulgrew B.Iterative frequency estimation by interpolation on Fourier coefficients ［J］.IEEE Transactions on Signal Processing，2005，53（4）：1237-1242.

［14］Provencher S.Estimation of complex single-tone parameters in the DFT domain ［J］.IEEE Transactions on Signal Processing，2010，58（7）：3879-3883.

［15］StoerJ，BulirschR.Introductiontonumericalanalysis［M］.New York：Springer，2002：43-49.

［16］韦岗，李道远，陈芳炯.基于逆积分方程的正弦波频率估计［J］.华南理工大学学报:自然科学版，2007，35（10）：147-151.Wei Gang，Li Dao-yuan，Chen Fang-jiong.Sinusoidal frequency estimation based on inversion integral equation［J］.Journal of South China University of Technology：Natural Science Edition，2007，35（10）：147-151.

［17］Aboutanios E.A modified dichotomous search frequency estimator ［J］.IEEE Signal Processing Letters，2004，11（2）：186-188.

［18］Quinn B G.Recent advances in rapid frequency estimation ［J］.Digital Signal Processing，2009，19（6）：942-948.