基于三次樣條插值的頻域頻偏估計算法研究

張少侃，梁中英，高心煒

（廣州海格通信集團股份有限公司，廣東 廣州 510663）

0 引 言

在目前的數字通信系統中，收發端載波頻偏的存在會使通信系統的性能急劇惡化。文獻[1-3]都對頻偏估計問題進行了深入的研究與討論，可見研究頻偏估計問題具有重要的意義。

目前，常見的數字通信系統中，數據輔助的頻偏估計算法由于其精度高、估計范圍廣，因而受到了更多的關注。數據輔助的頻偏估計算法主要分為兩大類，分別是基于時域分析的估計方法和基于頻域分析的估計方法。實際系統中常用的時域估計方法有Kay算法、Fitz算法、M&M算法以及L&R算法[4-7]。在時域算法中，綜合考慮估計范圍和估計精度，M&M算法是不錯的選擇，但其對信噪比比較敏感，即使在低信噪比下，估計性能也會惡化。在頻域算法中，基于快速傅里葉變換（Fast Fourier Transformation，FFT）變換對頻譜進行精確插值的頻偏估計算法比較常見[8]。此種算法頻偏估計范圍廣，對信噪比不敏感，常常受到眾多通信系統的青睞。但是由文獻[9]可知，由于FFT自身的柵欄效應，估計精度和有用數據個數以及FFT變換點數有關，因而對于估計精度要求比較高的系統常常難以滿足要求。二次拋物線頻域內插頻偏估計方法雖然使得估計精度有所提高，但此種方法因為只用到二次擬合,且沒有分區間段擬合，所以估計精度提高程度有限，仍然很難滿足很多系統的精度要求。

本文根據擬合理論，對FFT變換后的頻譜進行分段擬合。每段使用三次函數進行擬合，從而提高了擬合精度。另外，本文從基本代數理論出發，推導出了分段三次函數擬合的插值公式。仿真結果說明了三次樣條插值算法的估計誤差在低信噪比下仍能接近修正的克拉 - 美羅界（Modified Cramer-Rao Bound，MCRB），與二次插值算法相比，其估計誤差更接近修正的克拉-美羅界，且復雜度沒有提高，因而將會有更好的應用場景。

1 頻域頻偏估計方法

本文考慮物理傳輸信道模型為加性高斯白噪聲（Additive White Gaussian Noise，AWGN）信道，且接收端信號已經實現了位同步，則接收端收到的基帶離散信號可以表示為：

式中，T為符號周期時間；N為符號個數；sk為基帶信號的復數形式；Δf為接收基帶信號的頻率偏差；nk為加性復高斯白噪聲，其均值為0，方差為σ2。

FFT 是離散傅里葉變換（Discrete Fourier Transform，DFT）的一種快速算法，是工程實踐中可實現的方法。FFT頻域頻偏估計算法的基本思想是搜索接收端去調制信號頻譜圖的峰值。理論上去調制信號頻譜圖的峰值位置為零頻位置，當信號帶有頻偏時，峰值位置便為頻偏值的位置，因而峰值位置即對應于載波頻偏。FFT作為頻域頻偏估計算法，其不但頻偏估計范圍大，可以達到50%的符號速率，而且還能夠在很低的信噪比下工作。此外，由于使用快速傅里葉變換實現的復雜度低，因此在實際系統中得到了廣泛應用。DFT頻域估計算法和頻率估計分別可以表示為：

式中，Nfft為DFT變換點數；sk*是信號sk的復共軛形式；為所要估計的頻偏值。

2 三次樣條插值算法

三次樣條插值算法的基本思想是在利用基帶去調制信號頻譜圖峰值進行頻偏粗估計的基礎上，根據三次樣條插值公式，利用峰值及其左右兩條譜線的幅度進行細估。根據FFT得到頻譜能量的最大值A(km)及其對應的頻率點km，同時可以得到最大值A(km)左右相鄰的兩個值A(km-1)和A(km+1)，及其對應的頻率點km-1和km+1。根據這3點可以對構成的兩個區間段進行獨立的三次樣條插值曲線擬合，最后在兩個區間段中找到估計的最大值位置k以及下面給出三次樣條插值的定義及求解過程。

2.1 三次樣條插值定義

令x0=km-1，x1=km，x2=km+1，對應的y0=A(km-1)，y1=A(km)，y2=A(km+1)。樣條曲線S(x)是一個分段定義的光滑曲線，其分段函數曲線用三次多項式進行擬合，兩個分段函數可寫為：

式中，a0，a1，b0，b1，c0，c1，d0，d1代表8個未知系數。

2.2 三次樣條插值求解

根據前面的定義，步長h=x2-x1=x1-x0=1，將其帶入以下條件中。

（a）由Si(xi)=yi,i=0,1 可以推導出：

（b）由Si(xi+1)=yi+1,i=0,1可以推導出：

從式（5）到式（8）可知8個未知數有6個方程，還需要兩個方程才能求解出8個未知數。根據FFT變換的平滑性可知，曲線兩端沒有受到任何讓它們彎曲的力，即這樣便可以得到：

聯立式（5）到式（9）就可以求解出8個未知數，求解的結果為：

得到了8個未知系數后，便得到了兩段區間的光滑曲線函數表達式，根據三次函數的特性，可以求解出兩個區間段內的最大值點為：

對式（11）求解二次方程的根可得：

根據前面y0、y1、y2的定義可知，y1＞y0，y1＞y2，因此根號里面的分數值大于零。又因為三次曲線是建立在各段區間內的函數，而區間段長度h=1，所以對于式（13）取區間段內的根值為：

同時要求：

也即當y2＜y0時，在區間段[x0,x1)內，存在三次光滑曲線函數S0(x)的最大值，否則最大值不存在于區間段 [x0,x1)內。

對式（12）求解二次方程的根可得：

根據y1＞y0，y1＞y2可知，根號里面的分數值同樣大于零。同理，因為三次曲線是建立在各段區間內的函數，區間段長度h=1，所以對于式（16）取區間段內的根值為：

同時也還要求：

也即當y0＜y2時，在區間段[x1,x2)內，存在三次光滑曲線函數S1(x)的最大值，否則最大值不存在于區間段 [x1,x2)內。

2.3 三次樣條插值算法總結

綜上所述，下面給出三次樣條插值算法的頻域頻偏估計算法具體步驟。

步驟1：對接收到的基帶去調制信號進行快速傅里葉變換，找到頻域峰值坐標點(km,A(km))及其左右兩點 (km-1,A(km-1))與 (km+1,A(km+1))。

步驟2：判斷A(km-1)與A(km+1)的大小，如果A(km-1)＜A(km+1)，則轉入步驟3；如果A(km-1)＜A(km+1)，則轉入步驟4。

步驟3：此時計算估計出的最大值位置為：

轉入步驟5。

步驟4：此時計算估計出的最大值位置為：

轉入步驟5。

步驟5：計算頻偏估計值為：

3 仿真分析

對本文提出的基于三次樣條插值頻域頻偏估計算法進行MATLAB仿真分析。選取直接FFT頻域頻偏估計算法和二次插值頻域頻偏估計算法與本文提出的三次樣條插值算法做比較，對于FFT頻域頻偏估計算法：

對于二次插值頻域頻偏估計算法：

首先會對3種算法的頻偏估計范圍做仿真分析，研究它們在不同大小頻偏下的頻偏誤差均值，然后將頻偏估計結果與CRB進行比較。文獻[10]給出了信號載頻頻偏估計的MCRB，具體為：

式中，ES表示每個符號的能量；σ2表示復噪聲的功率。本文使用MCRB指標來對比算法的性能。

仿真設置調制方式為正交相移鍵控（Quadrature Phase Shift Keying，QPSK），假設接收有用的去調制符號個數N為64，信道為加性高斯白噪聲信道，初始相移由隨機函數產生，FFT的點數設置為256，蒙特卡洛實驗統計次數設置為1 000次。

3.1 頻偏估計范圍仿真分析

設置系統符號速率為4 ks/s，頻偏值的變化范圍為-2～+2 kHz，間隔為10 Hz，信噪比（Signal-to-Noise Ratio，SNR）=-10 dB，使用FFT算法、二次插值算法以及本文提出的三次樣條插值算法估計出的頻偏估計歸一化均方根誤差隨系統歸一化頻偏的變化曲線如圖1所示，歸一化步進為0.025。

圖1 不同算法的歸一化均方根誤差隨頻偏變化曲線

從上述仿真結果可以看出，3種算法均能在[-sr/2,sr/2]的頻偏范圍內工作，而且本文提出的三次樣條插值算法的歸一化均方根估計誤差在任何頻偏下都比FFT算法和二次插值算法的歸一化均方根估計誤差要小，有力地說明了本文給出的三次樣條插值算法在[-sr/2,sr/2]的頻偏范圍內均可以正常工作。

3.2 頻偏估計算法性能的仿真分析

設置系統符號速率為4 ks/s，頻偏值大小在[-sr/2,sr/2]的范圍內隨機產生，系統信噪比值的變化范圍為-10～+10 dB。FFT算法、二次插值算法以及本文給出的三次樣條插值頻偏估計算法的頻偏估計歸一化均方誤差和MCRB隨信噪比值變化的曲線如圖2所示。

圖2 不同算法的頻偏估計歸一化均方誤差與MCRB的比較曲線

從上述仿真結果可以看出，隨著信噪比的增加，各種算法的誤差曲線都在向MCRB曲線靠攏，但是在低信噪比下，三次樣條插值算法的誤差明顯更接近于MCRB曲線，且不會隨著信噪比的降低急劇增大。

4 結 論

針對FFT柵欄效應對頻域頻偏估計算法造成的精度問題，本文根據擬合理論以及基本代數理論建立了FFT變換后頻譜分段擬合模型，并通過詳細推導，求解出了三次樣條插值函數的最大點內插公式。最后，通過MATLAB仿真對比了FFT頻域頻偏估計算法、二次插值頻域頻偏估計算法以及三次樣條插值頻域頻偏估計算法的估計范圍和均方誤差，有效地說明了本文提出的三次樣條插值頻域頻偏估計算法的可行性與有效性。