專利名稱：一種基于多項式預(yù)測模型的多普勒頻率和相位估計方法
技術(shù)領(lǐng)域：
本發(fā)明屬于無線通信和衛(wèi)星導(dǎo)航接收機(jī)技術(shù)領(lǐng)域，具體涉及一種基于多項式預(yù)測
模型的多普勒頻率和相位估計方法。
背景技術(shù)：
由于多普勒效應(yīng)的存在，使得通信系統(tǒng)的收發(fā)載波可能存在較大的頻偏，導(dǎo)致接 收機(jī)無法解調(diào)出正確的發(fā)射數(shù)據(jù)[1]。特別是在以全球定位系統(tǒng)(GPS)為代表的擴(kuò)頻衛(wèi)星通 信中，載波跟蹤環(huán)路和碼相位跟蹤環(huán)路一起，構(gòu)成了接收機(jī)最重要的兩個組成部分[2'3]。
在過去的幾十年中，針對多普勒頻率的估計和跟蹤，文獻(xiàn)中提出了大量的方法。鎖 相環(huán)算法[2—5]由于具有結(jié)構(gòu)簡單、需要先驗知識少的特點而被廣泛采用，是一種對多普勒相 位的最大似然估計。另一類基于模型的、在最小均方誤差意義下最優(yōu)的多普勒頻率估計算 法包括擴(kuò)展卡爾曼濾波算法[6]，無氣味卡爾曼濾波算法[7—9]以及粒子濾波算法等等。在 多普勒頻率先驗信息已知的前提下，基于準(zhǔn)確模型算法的性能要明顯好于那些非基于模型 的算法。但是當(dāng)采用的模型不準(zhǔn)確時，算法的性能將惡化，甚至發(fā)散[11'12]。因此動態(tài)模型 的選擇，對這類基于模型的算法十分重要。 本發(fā)明首先提出一種全新的描述多普勒頻率和相位的動態(tài)模型——多項式預(yù)測 模型，在此模型的基礎(chǔ)上結(jié)合無氣味卡爾曼濾波算法，得到了一種新的估計多普勒頻率的 狀態(tài)空間方法——多項式預(yù)測_無氣味卡爾曼濾波方法。在狀態(tài)濾波的框架下，該算法還 可利用對其新息均值的檢測和判斷對模型進(jìn)行自確認(rèn)。無論收發(fā)機(jī)之間的相對運(yùn)動如何變 化，只要相對運(yùn)動的表達(dá)式滿足分段多項式形式，本發(fā)明提出的算法都可以有效的對多普 勒頻率和相位進(jìn)行估計。

發(fā)明內(nèi)容
本發(fā)明的目的在于提出一種應(yīng)用于GPS接收機(jī)的多普勒頻率和相位估計算法。
1.多普勒頻率和相位的多項式預(yù)測模型 接收機(jī)的多普勒頻偏是由于收發(fā)雙方的相對運(yùn)動引起的[1]，其表達(dá)式為
f = vA (1) 其中v為收機(jī)發(fā)機(jī)之間的相對速度，A為載波波長。通常，系統(tǒng)的載波波長是確 切已知的，這樣根據(jù)式(l)，接收機(jī)和發(fā)射機(jī)之間的相對運(yùn)動就完全決定了多普勒頻率。根 據(jù)牛頓運(yùn)動學(xué)公式，對于加速度恒定的運(yùn)動，其速度的表達(dá)式為
vt = v0+at (2) 其中vt為收發(fā)雙方在時間t時刻的相對速度，v。是初始時刻的相對速度，a為相 對運(yùn)動的加速度。 結(jié)合式(1)和式(2)，可以得到加速度恒定時多普勒頻率和相位的表達(dá)式為
f = v0/A+atA et= e 。+2 Ji tv0/A + Ji at2/A (3)
其中ft和e t分別是t時刻的多普勒頻率和相位，e 。是初始多普勒相位。由上式 可以看出，在加速度恒定的情況下，多普勒頻率和相位滿足多項式形式。當(dāng)式(3)被離散化 后，可以使用多項式預(yù)測濾波器[13'14]對式(3)建立狀態(tài)方程。
多項式預(yù)測濾波器的基本概念是對于一個離散化的L階多項式信號 = ，
(4)
/=0(其中p(l) (1=0，1，... ，L)為多項式的系數(shù)，n是離散的時刻)，可以用x(n)及其 前面(K-l)個時刻的值[x(n-K+l)，... ，x(n-l)]的線性組合來預(yù)測信號的將來值x (n+N)， 即 — + AO = Z ,傘_ " °
(5)
A=0 式(5)是一個以h(k) (k = 0， 1， . . . ， K-l)為系數(shù)的FIR濾波器，稱為多項式預(yù)測
濾波器?？紤]到實際系統(tǒng)中的信號都是混有噪聲的，若約束濾波器的噪聲增益S!l WWf最小。利用拉格朗日數(shù)乘法，結(jié)合式(4)和式(5)，使噪聲增益最小的多項式濾波器系數(shù)的最 優(yōu)解為[13]:
當(dāng)N = 1， L = 1時，
=- (6)
當(dāng)N = 1， L = 2時，

+ (_27 _ 36攀+ + 42A: +18
(7)
r-3K2 + 2^ 當(dāng)N和L為其它值時所對應(yīng)的最優(yōu)解可參見文獻(xiàn)[13]。值得注意到是，多項式預(yù) 測濾波器的階數(shù)需要滿足K^L+1。從式(6)和式(7)中可以看出，多項式預(yù)測濾波器式 的系數(shù)h (k)只與N、K和L有關(guān)，而與具體的描述信號的多項式系數(shù)p (1) (1 = 0， 1，. . . ， L) 無關(guān)，這意味著用多項式預(yù)測濾波器對多項式信號進(jìn)行預(yù)測時，不需要知道多項式信號本 身系數(shù)等先驗信息。基于式(5)的系統(tǒng)狀態(tài)方程可以表示為

xt — APPMxt—丄
其中AppM為
(8)A.
PPM
.Zz(O) 1 0
攀-1). 0
0 …1 0 h(k) (k = 0， 1， . . . ， K-l)為多項式預(yù)測濾波器的系數(shù)。 在本發(fā)明中，對多普勒頻率和相位同時進(jìn)行估計。多普勒頻率的狀態(tài)方程可以參 照式(8)的形式而直接得到；對于多普勒相位，基于恒加速度的假設(shè)，多普勒相位可以寫為 6 , = 6^ + w,—!7; + 0.5wTs2 = 6U + w,一i ; + 0.57;- w,_2], (9) 《+1.5^7;-0
5
其中9 t是t時刻的多普勒相位，Wt是t時刻的多普勒角頻率，Ts為采樣周期，w' 是角頻率變化率。 綜上，基于多項式預(yù)測模型的多普勒頻率和相位的狀態(tài)方程為

xt — ApPMXt—i，
其中Xt = [ e t
(10)
Wt Wt—J1. APPM的表達(dá)式為
A
ppm
1,57; -0.57; 畢)
1 0
(11)
+ n.
j方程)可以表示為
(12) h(l) (1 = 0， 1)是當(dāng)N = 1， L = 1， K = 2.時多項式預(yù)測濾波器的系數(shù)。
式(11)中的參數(shù)h(l)與式(3)中描述多形式的具體系數(shù)無關(guān)。這表明用式(5) 所描述的多項式預(yù)測濾波器式或式(8)所描述的狀態(tài)空間模型式來等效表示式(3)的離 散形式時，無需知道收發(fā)雙方確切的相對運(yùn)動的加速度、初始速度和初始多普勒相位等參 數(shù)信息，而不存在任何等效誤差。即在多項式階數(shù)確切已知的情況下，多項式預(yù)測模型式 (8)不需要引入描述模型不準(zhǔn)確性的過程噪聲，也不需要知道確切的描述多項式信號的多 項式系數(shù)，這與文獻(xiàn)[11]中描述的多項式模型是不同的。
2基于多項式預(yù)測模型的無氣味卡爾曼濾波算法
在加性高斯噪聲的信道中，假設(shè)接收到的信號(
—4 cos(《)-
z,=
L 4 sin(S,) 其中，~為接收到信號的幅度，9 t為t時刻的多普勒相位，nt為高斯白噪聲，其協(xié) 方差矩陣為R。 結(jié)合式(10) (11)和式(12)，可以得到一種新的狀態(tài)空間模型來對多普勒頻率和 相位進(jìn)行描述。由于觀測方程是非線性的，任何基于狀態(tài)空間模型的非線性濾波算法都可 以用來迭代求解，基于性能、運(yùn)算量和模型自確認(rèn)的綜合考慮，本發(fā)明選用無氣味卡爾曼濾 波算法進(jìn)行迭代求解。 由于多項式預(yù)測模型中不包含過程噪聲，因此在使用無氣味卡爾曼濾波算法進(jìn)行 估計時，設(shè)定過程噪聲的協(xié)方差矩陣為QPPM = 0。 在本文的模型中，由于假定觀測噪聲為加性高斯白噪聲，為了減少Sigma點的 數(shù)量，降低計算量，本發(fā)明采用了非擴(kuò)展形式的無氣味卡爾曼濾波算法[1°'15'16]，多項式預(yù) 測_無氣味卡爾曼濾波算法的具體步驟如下所示
1.初始化 文0 = E[x。]P。=五[(x。 — 、)(x。 一50)t ]其中x。是狀態(tài)向量的初始值，s。為初始狀態(tài)向
(13)
:的期望，P。是初始狀態(tài)向:
:的協(xié)方
差矩陣'


2.對t G {1，2，....①)，
(a)產(chǎn)生S i gma點Xw =! iM ±如+義)P,—!]
Xt—工為在t時的狀態(tài)向量，3^為狀態(tài)向量Xt—工的期望。
下面類似上面加橫線的記號均表示為期望。
(b).時間更新<formula>formula see original document page 7</formula>
(c). 二次產(chǎn)生Sigma點，
<formula>formula see original document page 7</formula>
(d)觀測更新
<formula>formula see original document page 7</formula> 在上面的步驟中，函數(shù)hO為觀測方程中的非線性函數(shù)，Zt為觀測向量，A為復(fù) 合尺度因子，nx是狀態(tài)向量的維數(shù)，Wi(m)為計算均值時對應(yīng)Sigma點的權(quán)重，w，為計算協(xié) 方差時對應(yīng)Sigma點的權(quán)重。
3多項式預(yù)測模型的自確認(rèn) 在第1節(jié)中，假設(shè)收發(fā)雙方的相對運(yùn)動滿足恒加速度形式。但在實際中，收發(fā)雙方 相對運(yùn)動的加速度可能是時變或時不變的，這意味著不可以直接把多項式預(yù)測模型應(yīng)用到 實際信號中。根據(jù)Weierstrass^逼近定理，在封閉區(qū)間內(nèi)的任意連續(xù)函數(shù)，都可以用一個 多項式以任意精度逼近。因此用多項式來表示收發(fā)機(jī)之間連續(xù)的相對運(yùn)動是合理的，特別 是當(dāng)選取一個合適的觀測時間窗時，用一個低階的多項式或分段多項式來描述收發(fā)機(jī)之間 的相對運(yùn)動可以得到令人滿意的精度，也就是說在實際應(yīng)用中，時變加速度的相對運(yùn)動也 可以用一個多項式或分段多項式來進(jìn)行近似。但是，對于分段多項式來說，各段多項式之間 存在不連續(xù)的點，在這些不連續(xù)點，加速度從一個值突變到另一個值，即多普勒頻率的變化 率從一個值突變到另一個值，信號從一種多項式形式進(jìn)入到另一種多項式形式。如果在這 些不連續(xù)點應(yīng)用多項式預(yù)測模型，將不能對信號進(jìn)行準(zhǔn)確的預(yù)測，即用不連續(xù)點前的采樣 點來預(yù)測不連續(xù)點后的信號是不準(zhǔn)確的。為了跟蹤時變多普勒頻率或者相位，必須想辦法 檢測出這些不連續(xù)點并加以修正。 對于一個最優(yōu)卡爾曼濾波器，進(jìn)入穩(wěn)態(tài)后，其新息序列是一個零均值的高斯白噪聲序列[18]，但是當(dāng)濾波器不是工作在穩(wěn)態(tài)時，新息序列的這一統(tǒng)計特性將不存在，也就是 說新息序列將不再是一個高斯白噪聲序列，其均值也不為O。根據(jù)分析，狀態(tài)方程(8)只能 以一定的階數(shù)來描述一個多項式信號，當(dāng)用式(8)描述一個分段多項式信號時，這些不連 續(xù)點將導(dǎo)致算法從穩(wěn)態(tài)轉(zhuǎn)移到過渡態(tài)，破壞新息序列的統(tǒng)計特性。文獻(xiàn)[18，19]中給出了 不同的用于檢測新息統(tǒng)計特性的方法。本發(fā)明中采用xs分布檢測法則來檢測新息序列的
均值是否為o，從而檢測出這些不連續(xù)點。多項式預(yù)測-無氣味卡爾曼濾波的新息可以表
示為(Zr^H)，新息的協(xié)方差矩陣為Pvv.因此，StTPvv—、T是一個有m個自由度的x 2變量，其 中m是新息的維數(shù)。在給定的置信度水平a下，StTPvv—5J與x。2的相對大小體現(xiàn)了新息 序列的均值是否為零。當(dāng)S,Tp^S^;^時，新息序列的均值不為零，意味著濾波算法沒有進(jìn) 入穩(wěn)態(tài)或者穩(wěn)態(tài)被打破；當(dāng)SfTPv/S,T^^,時，證明濾波算法已經(jīng)收斂到穩(wěn)態(tài)。
文獻(xiàn)[20]指出，當(dāng)所有關(guān)于模型的假設(shè)全都正確時，估計誤差的協(xié)方差矩陣會逐 漸收斂到一個常數(shù)陣，此時濾波算法進(jìn)入穩(wěn)態(tài)。 一旦穩(wěn)態(tài)被破壞，使用狀態(tài)方程得到的狀態(tài) 預(yù)測不再準(zhǔn)確，大量"準(zhǔn)確"的信息包含在觀測序列中，即包含在新息中。為了使算法能夠快 速再次收斂而不至于發(fā)散，新息增益矩陣需要被增大。在狀態(tài)方程和觀測噪聲方差確定的 情況下，新息增益隨著過程噪聲方差的增加而增加[20]。另一方面，過程噪聲的方差描述了 模型的不準(zhǔn)確性，當(dāng)模型不準(zhǔn)確時，過程噪聲的方差也應(yīng)該被增大。因此在這些不連續(xù)點， 當(dāng)穩(wěn)態(tài)被打破時，即S/Pvv"S,T:^J，狀態(tài)方程不再準(zhǔn)確，此時令Qppm = 13 I o n2，(其中o n2是 觀測噪聲的方差，13為一個大于l的正數(shù))，可以加速濾波算法的收斂，使其快速進(jìn)入下一 個穩(wěn)態(tài)。而當(dāng)濾波算法再次進(jìn)入穩(wěn)態(tài)后，即S,TPw"S/^A重新令Qppm = O，進(jìn)行正常的迭代
濾波。因此完整的多項式預(yù)測-無氣味卡爾曼濾波算法應(yīng)該包括自確認(rèn)的步驟來處理分段 恒加速運(yùn)動中的不連續(xù)點。在觀測更新之后，計算Sz-(Zr^H)并在一定的置信度水平a下
比較StTPvv—、t與x 2的關(guān)系。gS,TPw 、T〉;^時，令Op = 13 I 。 n2 ;當(dāng)S7Pw"S/時，令QPPM =0. 需要注意的是，當(dāng)相鄰兩個不連續(xù)點之間的時間間隔小于算法的收斂時間時，本 發(fā)明提出的方法將一直處于過渡態(tài)而無法進(jìn)入穩(wěn)態(tài)。此時算法退化為狀態(tài)模型不準(zhǔn)確且包 含較大過程噪聲的無氣味卡爾曼濾波算法。狀態(tài)估計結(jié)果主要依賴于新息對狀態(tài)預(yù)測結(jié)果 的更新，算法不會發(fā)散，但也無法達(dá)到穩(wěn)態(tài)時的最優(yōu)。
4.GPS調(diào)制信號的處理 在實際的GPS信號中，導(dǎo)航數(shù)據(jù)通過BPSK方式調(diào)制在載波上，這樣在載波上會不
定時的發(fā)生±180°相位翻轉(zhuǎn)。因此需要對數(shù)據(jù)進(jìn)行一定的處理，克服相位翻轉(zhuǎn)對頻率估計
的影響。在忽略噪聲的前提下，設(shè)接收到的1、Q兩路的信號滿足
/ = D0), /I cos(外))=D( ;M cos( w" p) (25)
g = D(f) 乂 sin(柳)=￡>0).爿sinOf + p) (26)
其中，A是載波的幅度，假設(shè)其為一常數(shù)，D(t)為調(diào)制信息，取值為±1， w為多普 勒頻偏。對于本方法提出的算法，為了克服±180°相位翻轉(zhuǎn)的影響，可以將觀測轉(zhuǎn)換為
D2 (t) A2cos [2 9 (t) ] = D2 (t) A2 [cos2 (9 (t)) -sin2 (9 (t)) ] = I2-Q2 (27)
D2(t)A2sin[2 9 (t) ] = 2D2 (t) A2sin ( 9 (t)) cos( 9 (t)) = 2IQ (28)
上式中D2(t) = l，這樣就避免了調(diào)制數(shù)據(jù)對載波相位的影響。觀測方程可以表示
為
—乂cos(2《)1 。o、 需要注意的是，當(dāng)采用式(29)作為觀測方程時，直接估計的結(jié)果為實際頻偏的2 倍，需要進(jìn)行后處理。同時為了避免信號幅度對估計結(jié)果的影響，需要對1、Q兩路信號進(jìn)行 了歸一化處理，即I = I/(I2+Q2) ， Q = Q/(I2+Q2)。
技術(shù)效果 本發(fā)明為一種基于多項式預(yù)測模型、應(yīng)用于GPS接收機(jī)中的估計多普勒頻率和相 位的算法法。無論收發(fā)機(jī)之間的相對運(yùn)動如何變化，只要相對運(yùn)動的表達(dá)式滿足分段多項 式形式，本發(fā)明提出的算法都可以有效的對多普勒頻率和相位進(jìn)行估計。


圖1仿真多普勒頻率信號隨時間變化的軌跡。 圖2本發(fā)明算法與所有對比算法在仿真信號情況下的性能比較，其中圖2a為在 17dB時，100次蒙特卡洛仿真的平均頻率估計誤差；圖2b為不同信噪比下，各種算法平方根 均方誤差的比較。 圖3為高動態(tài)多普勒頻率信號的軌跡。其中，3a為多普勒頻率信號，圖3b為頻率 一介變化率的，圖3c為頻率二介變化率的。 圖4本發(fā)明算法與對比算法在高動態(tài)多普勒頻率信號情況下的性能比較，其中圖 4a為17dB時，100次蒙特卡洛仿真的平均頻率估計誤差；圖4b為不同信噪比下，各種算法 平方根均方誤差的比較。 圖5本發(fā)明方法與對比算法對實際GPS信號的頻率估計實驗結(jié)果，其中圖5a為整 個實驗過程的頻率跟蹤結(jié)果；圖5b為在實驗初始階段頻率跟蹤結(jié)果。
具體實施例方式1.根據(jù)式(10) (11)建立多普勒頻率和相位的多項式預(yù)測模型； 2.若載波上無調(diào)制信息，則采用式(12)作為觀測方程；若載波上存在調(diào)制信息，
則利用式(27)、 (28)對接收信號進(jìn)行預(yù)處理，采用式(29)作為觀測方程。 3.令QppM = O，設(shè)定一定的置信度水平a和常數(shù)|3 ; 4.使用式(14-24)所示的無氣味卡爾曼濾波算法進(jìn)行迭代濾波，在每次迭代的最 后，計算新息的方差StTPvv—、t，如果大于閾值x a2，則令Op = 13 I o n2，如果小于閾值，則令 QppM — 0。 5.若無調(diào)制信息，則估計出的相位和頻率就是實際的多普勒相位與頻率，如果有 調(diào)制信息，則估計出的多普勒頻率和相位是實際相位和頻率的2倍，需要進(jìn)行后處理。
仿真結(jié)果
1.仿真條件 為了評估本發(fā)明所提出方法的性能，將本發(fā)明提出的方法與鎖相環(huán)算法、基于傳
9<formula>formula see original document page 10</formula>
統(tǒng)白噪聲激勵模型的擴(kuò)展卡爾曼濾波算法、無氣味卡爾曼濾波算法和粒子濾波算法進(jìn)行比
較。在下面的仿真中，為了表示方便，假設(shè)式(12)中接收信號的幅度為一個常數(shù)，不失一般 性，令觀測方程中At = 1。 對于本發(fā)明提出的算法，基于大量的數(shù)值仿真結(jié)果，在做新息均值檢測時，設(shè)定置 信度水平a為0.25，常數(shù)13 = 10。 對于鎖相環(huán)算法，在仿真和實驗中均采用了常用的二階鎖相環(huán)結(jié)構(gòu)[21]，鎖相環(huán)的 參數(shù)設(shè)置與文獻(xiàn)[5]中一致，噪聲帶寬Bn = 50Hz，阻尼比為0. 707. 在機(jī)動目標(biāo)運(yùn)動模型中，常用的白噪聲激勵模型包括白噪聲加速度(WNA)[11]模型
和白噪聲加加速度(WNJ)""模型。根據(jù)式(3)可知，機(jī)動目標(biāo)的速度與多普勒頻率存在著
線性映射關(guān)系，因此與WNA和WNJ所對應(yīng)的多普勒頻率的動態(tài)模型稱為白噪聲多普勒頻率
變化率模型和白噪聲多普勒頻率二階變化率模型。 白噪聲多普勒頻率變化率模型的狀態(tài)方程可以表示為 x(t) = Ax(t-l)+q(t-l) (30) 其中狀態(tài)向量x(t) = [ 9 a)W(t)]T，q(t-l)的協(xié)方差矩陣為Q:

，TS為系統(tǒng)的采樣時間，o v2為白噪聲多普勒頻率變化率的方差。
白噪聲多普勒頻率二階變化率的狀態(tài)方程可以表示為 x(t) = Ax(t-l)+q(t-l) (27)
其中狀態(tài)向量x(t) = [9 (t) w(t) w'a)]1, q(t-l)的協(xié)方差矩陣為
A =
o/是白噪聲多普勒頻率二階變化率的方
差。在下面的第一組和第二組仿真中分別用到了上述兩種模型。 在常溫下，熱噪聲的功率譜密度為_204dB W/Hz，對于GPS Ll信號，信號的帶寬 為2M Hz^，這樣熱噪聲功率約為-141dB W。對于一般的GPS接收機(jī)，要求的最低靈敏度 為-160dB wW，在最低靈敏度下，對應(yīng)的信噪比大約為-19dB，載噪比為44dB-Hz。由于GPS 信號采用擴(kuò)頻調(diào)制，擴(kuò)頻碼的碼片速率為1. 023M Hz，相關(guān)積分累加時間為lms，這樣就存在 著大約33dB的積分增益。因此在下面的各項仿真結(jié)果中，均給出了信噪比從14dB到20dB 的頻率估計結(jié)果，以符合實際系統(tǒng)的要求，而采樣時間Ts均設(shè)定為lms。
1. l多普勒頻率仿真信號 為了驗證本發(fā)明算法在信號滿足分段多項式形式和多普勒頻率滿足正弦連續(xù)變
化時的性能，給出下列一段多普勒頻率仿真信號
x2
V"2 5"
10<formula>formula see original document page 11</formula>
(30)
<formula>formula see original document page 11</formula> 這段信號中包含了多普勒頻率變化率為常數(shù)、隨機(jī)變量、連續(xù)變化以及發(fā)生突變 時的情況。式(30)中頻率變化率f' =25Hz/s， v是一個滿足正態(tài)分布的隨機(jī)變量，其方 差o2 = (25Hz/s)2。對應(yīng)多普勒頻率的軌跡如圖1所示。 在本組仿真中，由于多普勒頻率變化率隨時間而變化，因此基于白噪聲激勵模型 的對比算法均采用白噪聲多普勒頻率變化率模型，多普勒頻率變化率的方差ov2= o2 = (25Hz/s)2。 1. 2高動態(tài)多普勒頻率信號 為了驗證本發(fā)明算法在高動態(tài)模型下的性能，仿真使用了文獻(xiàn)[5]中的給出的高 動態(tài)模型，初始速度為_40111/8，在跟蹤過程中，有一恒定的加速度_25g(g為重力加速度)， 中間還有兩段各持續(xù)0. 5s，大小為100g/s2的加加速度。根據(jù)式(3)，在GPS高動態(tài)運(yùn)動模 型下，對應(yīng)的多普勒頻率高動態(tài)模型如圖3所示。 由于在本組仿真中，絕大部分時間里多普勒頻率變化率都為一個不等于零的常 數(shù)，因此在該仿真中采用白噪聲多普勒頻率二階變化率模型，多普勒頻率二階變化率的標(biāo) 準(zhǔn)差取為實際最大多普勒頻率二階變化率的三分之一，即3 o j = 4900Hz/s2， o " 1633Hz/ 1.3實際GPS信號 為了驗證本發(fā)明算法在實際系統(tǒng)中的有效性，采用實際接收到的GPS數(shù)據(jù)來對算 法進(jìn)行驗證。該組仿真中所用的實際GPS信號是文獻(xiàn)[4]隨書附帶光盤中的數(shù)據(jù)。該數(shù)據(jù)由 NordNav公司的模擬前端R30對GPS信號進(jìn)行接收采樣，采樣地點位于意大利的都靈。R30 模擬前端的采樣頻率為16. 3676M Hz，中頻頻率為4. 1304MHz，量化精度為4bit。所用數(shù)據(jù) 先經(jīng)過文獻(xiàn)[4]所提供軟件接收機(jī)的捕獲處理，并經(jīng)過粗頻偏估計和補(bǔ)償，使得殘余多普 勒頻偏的絕對值小于250Hz，并按照式(27) (28)對數(shù)據(jù)進(jìn)行預(yù)處理，消除調(diào)制數(shù)據(jù)的影響。 對經(jīng)過預(yù)處理后的數(shù)據(jù)，分別采用本發(fā)明中的方法和另外四種對比算法對多普勒頻率進(jìn)行 估計和跟蹤。在實驗中。我們截取了 2.8s的實際GPS數(shù)據(jù)，相干累積積分時間為lms?；?于白噪聲激勵模型的算法均采用白噪聲多普勒頻率變化率模型，多普勒頻率變化率的方差 o v2 = (25Hz/s)2。
2.實驗結(jié)果 在第一組仿真中，在多普勒仿真信號的各個階段，本發(fā)明提出的算法估計性能都 要好于其它對比算法，仿真結(jié)果如圖2所示。在t = 600時，多普勒頻率變化率發(fā)生突變，對 于本發(fā)明提出的模型，此時分段多項式出現(xiàn)不連續(xù)點。按照上一節(jié)的分析，此時本發(fā)明的方 法會檢測到不連續(xù)點，并在一段時間內(nèi)增大過程噪聲方差，直到算法重新收斂。因此當(dāng)t = 600時，本發(fā)明提出方算法有一個估計誤差的小的峰值，其后迅速下降，算法再次收斂，從而 說明了本方面方法中自確認(rèn)部分的有效性。當(dāng)多普勒頻率滿足正弦變化時，由于不是嚴(yán)格 滿足多項式形式，本發(fā)明算法的性能有所下降，但還是明顯好于其它對比算法。
高動態(tài)多普勒頻率信號情況下的仿真結(jié)果如圖4所示。從途中可以看出，本發(fā)明提出的算法在多普勒頻率二階變化率為0時，明顯好于其它對比算法。而當(dāng)多普勒頻率二階變化率的絕對值為一個很大的值時，由于多普勒頻率變化率不斷發(fā)生變化，用本發(fā)明提出的基于分段恒加速度假設(shè)的模型將不能準(zhǔn)確的描述多普勒頻率變化率的快速變化情況。因此在區(qū)間t G [2000,2500]和t G [4500， 5000]，本發(fā)明提出的方法一直處于過渡態(tài)，過程噪聲的值較大，算法退化為模型不準(zhǔn)確時無氣味卡爾曼濾波算法，其估計性能要比采用白噪聲激勵模型的方法要差，但是在穩(wěn)態(tài)時，本發(fā)明提出的算法性能要明顯好于其它對比算法。 實際信號的頻率跟蹤結(jié)果如圖5所示。從圖中可以看出，各種算法對于實際載波
信號的頻率估計都可以較快收斂到一個比較穩(wěn)定的值，而不會發(fā)生發(fā)散，并且本發(fā)明提出
算法的估計結(jié)果的方差明顯要小于其它對比算法，從而驗證了本文算法的有效性，可以在
實際的GPS接收機(jī)中得到應(yīng)用。 參考文獻(xiàn) [l]Proakis J G， Digital Communications[M]，3rd ed. New York :McGrawHill，1995. [2]E D Kaplan and C. J. Hegatry. Understanding GPS Principles andApplications, ed. Norwood, MA :Artech House, 2006. [3]James Bao_Yen and TSUI. Fundamentals of Global Positioning SystemReceivers-Asoftware Approach,2nd ed. John Wiley & Sons, Inc. Hoboken， New Jersey,2005. [4]Kai Borre，Dennis M. Akos，Nicolaj Bertelsen et al. A Software-definedGPS andGalileo Receiver. Boston :Birkhauser，2006. [5]V A Vilnrotter，S Hinedi，R Kumar. Frequency Estimation Techniques forHigh DynamicTrajectories. IEEE Trans, on Aerospace and Electronic Systems. 25(4):1989. [6]Steven M. Kay.Fundamentals of Statistical Signal Processing VolumeI :EstimationTheory. New Jersey :Prentice Hall，1993. [7] Simon J. Julier， Jeffrey K. Uhlma皿.A New Extension of the KalmanFilterto NonlinearSystems，Proc. Of AeroSense :The 11th International Symposium onAerospace/Defence Sensing, Simulation and Controls, Orlando, Florida, Vol. MultiSensor Fusion, Tracking and ResourceManagement II，1997. [8]Rudolph van der Merwe，Arnaud doucet，nando de Freitas and Eric Wan.The皿scentedparticle filter, Aug. 2000，Oregon Graduate Institute
[9]Julier， Simon J， Jeffrey K. Uhlman. Unscented filtering and nonlinearestimation. Proceedings of the IEEE，vol.92(3) :401—422，2004. [10]M. Sanjeev Arulampalam， Simon Maskell， Neil Gordon and Tim Clap.A tutorial o即article filters for online nonlinear non_Gaussian Bayesiantracking. IEEE Trans. On SignalProcessing，vol. 50(2) :174-188，2002.
[ll]X. Rong Li and Vesselin P. Jilkov， Survey of Maneuvering TargetTracking. Part I :Dynamic Models. IEEE Trans, on Aerospace and Electronic Systems,2003,39(4) :1333-1364. [12]Bar_Shalom， Y. ， Li， X. R. ， and Kirubarajan， T. (2001)Estimation withAppliations toTracking and Navigation:Theory，Algorithm,and Software. New York:Wiley，2001. [13]P. Heinonen and Y. Neuvo， FIR-Median Hybrid Filters with PredictiveFIRSubstructures [J]， IEEE Trans.on Acousitics， Speech and Signal Processing,1998,36(6) :892-899. [14] Sami Valiviita， Seppo J. Ova ska and Olli Vainio. PolynomialPredictive Filtering inControl Instrumentation :A Review. IEEE.Trans. OnIndustrial Electronics, vol.46 (5)，1999. [15]Jouni HartiKainen， Simo Sarkka.Optimal filtering with Kalmanfilters and smoothers-a Manual for Matlab toolbox EKU/UKF. [16]Wu， Y. ， Hu， D， Wu M and Hu X. Unscented Kalmen Filtering for additivenoise case :augmented versus nonaugmented. IEEE Signal Processing Letters 12(5):357-360,2005. [17]Eric W Weisstein. 〃 Weierstrass Approximation Theorem. 〃 FromMathWorld—AWolfram Web Resource. http://mathworld. Wolfram, com/WeierstrassApproximationTheorem. html [18]R. K. Melrra， On the identification of variance and adaptive Kalmanfiltering. IEEETrans. OnAutomatic Control, Vol 15(2) :175-184,1970.
[19]E. Parzen. An approach to time series analysis. Ann. Math， Statist,Vol. 32 (5) :951-989,1961. [20]Fitzgerald R. J. ， Divergence of the Kalman Filter, IEEE. Trans, onAutomation Control, 1971，16(6) :736-747. [21]Roland E. Best. Phase-Locked Loop Design,Simulation and Application,4th edition. Columbus, OH :McGraw-Hill， 1999.
1權(quán)利要求
一種基于多項式預(yù)測模型的多普勒頻率和相位估計方法，其特征在于采用描述多普勒頻率的動態(tài)模型——多項式預(yù)測模型，在恒定多普勒頻率變化率的假設(shè)條件下，基于多項式預(yù)測模型的多普勒頻率和相位的狀態(tài)方程表示為xt＝APPMxt-1，其中xt＝[θt wt wt-1]T，θt為t時刻的多普勒相位，wt為t時刻的多普勒角頻率，APPM的表達(dá)式為h(l)，l＝0，1，是當(dāng)N＝1，L＝1，K＝2，時多項式預(yù)測濾波器的系數(shù)，Ts為采樣周期，L為多項式預(yù)測濾波器的多項式信號的階數(shù)。F2010100230630C00011.tif
1. 一種基于多項式預(yù)測模型的多普勒頻率和相位估計方法，其特征在于采用描述多普 勒頻率的動態(tài)模型——多項式預(yù)測模型，在恒定多普勒頻率變化率的假設(shè)條件下，基于多項式預(yù)測模型的多普勒頻率和相位的狀態(tài)方程表示為 Xt — APpMxt_i，其中Xt = [ e t Wt Wt—jt， e t為t時刻的多普勒相位，wt為t時刻的多普勒角頻率，AppM 的表達(dá)式為—l 1.57; -0.5r/0 畢)Zz(l) 0 1 0h(l)，l = O，l，是當(dāng)N二 l，L = 1，K = 2，時多項式預(yù)測濾波器的系數(shù)，Ts為采樣周期， L為多項式預(yù)測濾波器的多項式信號的階數(shù)。
2. 根據(jù)權(quán)利要求1所述的多普勒頻率和相位估計方法，其特征在于在加性高斯噪聲的 信道中，假設(shè)觀測方程表示為—4 cos(W 4 sin(61,)其中，At為接收到信號的幅度，e t為t時刻的多普勒相位，nt為高斯白噪聲，其協(xié)方差 矩陣為R ;將多項式預(yù)測模型與無氣味卡爾曼濾波算法相結(jié)合進(jìn)行迭代求解，具體步驟如 下1)初始化:[50 = E[x0]P0 = ￡[(x0 -Io)(x。 -、)T] (13)其中x。是狀態(tài)向量的初始值J。為初始狀態(tài)向量的期望，P。是初始狀態(tài)向量的協(xié)方差矩(12)陣2)對t G {1，2，(a)產(chǎn)生Sigma點XM ± V(似+義)P,-1 ]xt—i為在t時的狀態(tài)向量，^為狀態(tài)向量Xt—工的期望 下面類似上面加橫線的記號均表示為期望；時間更新—ApPMxt_i(14)2wx(15)2nxp'k-i = ZXc)[x'，,i,-i - Vi][x"l'-1 - Vi]t +Qf*=0(c) 二次產(chǎn)生Sigma點，yt|t-1(d)=h(x' t|t—"t)2 xy，1 = 2>;m)y,,,i,—1l測更新pv,巾-i - Vi][yM'-1UT+R,'=02(16) (17)(18)(19)(20)2ra:(21)K, = P P_1& = Vi+K'(z'—1) P' =IVi - K,pv,v,K,在上面的步驟中，函數(shù)h(—為)(22) (23) (24)J方程中的非線性函數(shù)，zt為;l向量，a為復(fù)合尺度因子，nx是狀態(tài)向量的維數(shù)，w，為計算均值時對應(yīng)Sigma點的權(quán)重，w，為計算協(xié)方差 時對應(yīng)Sigma點的權(quán)重。
3. 根據(jù)權(quán)利要求2所述的多普勒頻率和相位估計方法，其特征在于在狀態(tài)濾波的框架 下，在每次迭代最后，利用對新息均值的檢測和判斷進(jìn)行模型的自確認(rèn)；即設(shè)無氣味卡爾曼 濾波器的過程噪聲協(xié)方差矩陣為QppM，新息為St，新息的協(xié)方差矩陣為Pvv，在給定的置信度 水平a下，當(dāng)StTPvv-、t > x a2時，QPPM = |3 I 0 n2，當(dāng)StTPvv-、t《x a2時，QPPM = 0。
4. 根據(jù)權(quán)利要求2所述的多普勒頻率和相位估計方法，其特征在于當(dāng)載波上存在周期 信息，則對信號進(jìn)行預(yù)處理后，使用如下觀測方程Jcos(2(9,)-爿sin(2(9,)
全文摘要
本發(fā)明屬于無線通信和衛(wèi)星導(dǎo)航接收機(jī)技術(shù)領(lǐng)域，具體為一種基于多項式預(yù)測模型的多普勒頻率和相位估計方法。本發(fā)明首先提出一種新的描述多普勒頻率和相位的動態(tài)模型——多項式預(yù)測模型，并在此模型基礎(chǔ)上結(jié)合無氣味卡爾曼濾波器，提出一種新的自確認(rèn)濾波算法用以估計多普勒頻率和相位。無論收發(fā)機(jī)之間的相對運(yùn)動如何變化，只要相對運(yùn)動的表達(dá)式滿足分段多項式形式，本發(fā)明提出的方法都可以有效的對多普勒頻率和相位進(jìn)行估計。數(shù)值仿真結(jié)果和真實GPS接收機(jī)實驗數(shù)據(jù)均表明，在GPS接收機(jī)可正常工作的信噪比條件下，本發(fā)明提出的方法優(yōu)于大部分其他文獻(xiàn)中已報道的多普勒頻率估計算法。
文檔編號H04L27/26GK101753513SQ20101002306
公開日2010年6月23日 申請日期2010年1月21日 優(yōu)先權(quán)日2010年1月21日
發(fā)明者尹建君, 張建秋, 趙晉 申請人:復(fù)旦大學(xué)