去中心化時差頻差直接定位方法

朱穎童， 董春曦， 董陽陽， 許錦， 趙國慶

西安電子科技大學 電子信息攻防對抗與仿真技術教育部重點實驗室， 西安 710071

去中心化時差頻差直接定位方法

朱穎童*， 董春曦， 董陽陽， 許錦， 趙國慶

西安電子科技大學 電子信息攻防對抗與仿真技術教育部重點實驗室， 西安 710071

針對原始利用時差頻差的直接定位(DPD)方法存在數據傳輸量和計算量大的瓶頸，提出了兩種去中心化直接定位方法。第1種方法采用去中心化配對方案，只將各觀測站截獲信號在站間進行一次傳輸，將數據傳輸和計算分散到各觀測站間并行計算互模糊函數(CAF)，構造僅滿足滿秩條件的互模糊矩陣(CAM)。第2種方法根據推導的任意互模糊函數間關系公式，采用歸約方式去中心化的在各觀測站并行計算余下互模糊函數，補全互模糊矩陣。兩種方法都降低了直接定位數據傳輸量，提高了計算效率。性能分析和仿真實驗表明本文兩種方法精度性能優于兩步定位方法，在低信噪比時兩種方法都可達到比較理想的精度性能，在高信噪比時第2種方法與原始直接定位方法的精度性能相當。

直接定位； 無源定位； 互模糊函數； 去中心化； 觀測站配對； 克拉美羅下界

無源定位具有作用距離遠，安全隱蔽性能好等優勢，是電子對抗的一個重要的研究方向[1]。無源定位已廣泛應用于飛機、衛星、無線傳感器網絡等平臺。近年來將無源定位應用于低成本無人機蜂群(LOw-Cost UAV Swarm Technology，LOCUST)這類自組協同受尺寸、重量、功率、成本(Size Weight and Power Cost，SWaP-C)[2]約束的分布式平臺，對定位算法提出了更苛刻的要求。無源定位一般采用兩步定位方法，首先利用截獲信號估計到達角(Direction Of Arrival，DOA)、到達時間(Time Of Arrival，TOA)、時差(Time Difference Of Arrival，TDOA)、頻差(Frequency Difference Of Arrival，FDOA)和信號強度(Received Signal Strength，RSS)等參數，再利用這些參數估計輻射源位置。

Wax和Kailath[3]首先提出了直接利用多個陣列的信號采樣，計算協方差矩陣實現窄帶輻射源數目及位置的一步估計，并討論了位置與頻率聯合估計以及運動輻射源的情況。Weiss等[4-7]將直接利用信號實現輻射源位置一步估計的方法命名為直接定位(Direct Position Determination，DPD)方法。這類方法避免了兩步方法所需的信號分選和參數配對過程，且可在低信噪比下得到精度更高的定位結果。Oispuu和Nickel[8]首先提出運動單站天線陣列對多個間歇輻射的輻射源目標的檢測和直接定位方法。郭福成等[9-10]提出利用運動單站干涉儀接收的信號數據估計輻射源位置的直接定位方法，采用一組干涉儀便可對多個頻段輻射源進行定位，且在低信噪比下優于多通道干涉儀僅測角和長基線相位差變化率定位方法。王鼎等[11]也采用運動單站天線陣，依據相位調制信號具有的恒包絡特性，根據最大似然準則實現輻射源的直接定位，明顯提高了位置估計的精度。Bosse等[12]提出Global MUSIC方法實現多站天線陣對多個輻射源目標的直接定位。黃志英等[13]在實現多站天線陣對多個輻射源直接定位時考慮陣列的非一致噪聲，通過交替迭代的方式對噪聲協方差矩陣以及目標位置進行聯合估計。

Weiss等在文獻[6-7]中分別提出了多個觀測站對單個固定輻射源的直接定位方法，每個觀測站只有一個接收機截獲一路信號，基于時差或頻差或兩者聯合,分別適用于窄帶和寬帶的未知或先驗已知信號。Kay等[14-15]利用多個觀測站實現輻射源的低截獲(Low Probability of Intercept，LPI)信號檢測和直接定位，針對采樣信號存在的信號量化效應影響定位性能，采用了連續時間信號模型，適用于窄帶、寬度、低通和高通的輻射時間未知的未知或已知信號。直接定位方法由于利用信號進行定位，因此存在著數據傳輸量大、計算量大等瓶頸。

為了克服這些瓶頸，Fowler等在文獻[2]和文獻[16-17]中分別提出了數據壓縮和分布式計算的方法，并分析了其對直接定位精度性能的影響。王云龍和吳瑛[18]利用矩陣轉置后的非零特征值相同這一性質，降低代價函數中矩陣的維數，減小其特征值分解的計算量。夏威等[19-20]提出自適應分布式計算實現輻射源直接定位，將各觀測站截獲采樣的信號傳輸給鄰近觀測站而非單一觀測站或者融合中心(Fusion Center)，顯著降低了數據通信量且將計算量分散到各觀測站。另外采用最速下降的迭代方式，不需要遍歷計算所有格點，降低了估計輻射源位置所需的總計算量。

相關基于時差或頻差或兩者聯合的直接定位，需將L個觀測站兩兩配對成共L(L-1)/2組配對，將各站截獲信號兩兩配對逐一計算廣義互相關(Generalized Cross Correlation，GCC )[21]或者互模糊函數(Cross Ambiguity Function，CAF)[22]，每個網格點再構造包含輻射源位置信息的Hermite矩陣的L(L-1)/2個矩陣元素。Hermite矩陣即為互模糊矩陣(Cross Ambiguity Matrix，CAM)。而兩步定位法只利用L-1組配對，提取L-1組時差和頻差用于定位解算。另外在兩步定位方法中，由于中間變量包含有輻射源與某一觀測站之間的斜距，需要設定主觀測站和副觀測站，所用參數一般為主副觀測站間的時差或頻差。定義主副觀測站間的時差和頻差為中心化時差和頻差。

由于各觀測站截獲的是相同輻射源同一時刻的信號，GCC或CAF間滿足相應關系。因此，提取L-1組時差頻差，只計算L-1個GCC或CAF，其他(L-1)(L-2)/2個函數不需計算GCC或CAF，可利用函數間關系公式產生，以補全互模糊矩陣，極大降低了觀測站間數據傳輸量和計算量。

綜上，本文嘗試提出兩種方法構造互模糊矩陣，第1種方法不需全部計算L(L-1)/2個矩陣元素，只需要滿足矩陣滿秩的條件，最少只需要計算L-1個矩陣元素，降低計算量。與文獻[2，19-20]類似，不設定主觀測站或者融合中心，采用去中心化的觀測站配對方案，觀測站只向鄰近的觀測站轉發截獲信號，極大降低了定位系統的數據傳輸量。采用分布式的計算方式，將原來相關運算都在同一個觀測站的方式改為分散到各個觀測站分別計算，縮減了定位計算所消耗的總時間。第2種方法在第1種方法計算得到的L-1個GCC或CAF基礎之上，利用函數之間的關系計算余下(L-1)(L-2)/2個GCC或CAF，在保證構造出所有矩陣元素的情況下，極大降低了數據傳輸量和計算量。

由于時差頻差聯合的直接定位可簡化得到時差或頻差的直接定位，因此本文討論時差頻差聯合的直接定位。本文內容組織如下：第1節給出模型并總結相關直接定位方法，第2節給出兩種互模糊矩陣構造方法，第3節推導各定位方法克拉美羅下界(Cramer-Rao Lower Bound，CRLB)，第4節進行仿真實驗，第5節給出結論。

1 問題模型

假設輻射源目標位置為p。L個觀測站的位置和速度分別為pl和vl，l=1,2,…,L。各觀測站K個時間段T內截獲的輻射源信號表示為

rl,k(t)=bl,ksk(t-tl,k)ej2πfl,kt+wl,k(t)

0≤t≤T

(1)

式中：k=1,2,…,K；bl,k為信號輻射到達觀測站l的傳播損耗；sk(t)為信號波形；wl,k(t)為零均值加性高斯白噪聲；T為觀測時間；tl,k和fl,k分別為第k時間段內信號到達觀測站l的到達時間和多普勒頻移，其表達式分別為

(2)

(3)

信號到達觀測站l和m的時差和頻差分別為

(4)

flm,k=fl,k-fm,k=

(5)

為了討論的簡便，假定在觀測時間內滿足τlm,k?T，即定位過程中時差和頻差保持不變。

文獻[16-17]指出對截獲的信號進行采樣和壓縮傳輸會產生量化誤差，文獻[15]指出觀測站間的同步誤差影響估計結果，文獻[23]指出運動觀測站間截獲信號的轉發傳播誤差影響估計結果。

文獻[2,6,14,19-20]對各觀測站截獲信號時域直接采樣，采樣率限制了時差的量化精度。文獻[16-17]將信號進行離散傅里葉變換、小波包變換和短時傅里葉變換后再量化壓縮，壓縮時在時差頻差精度，計算量與傳輸量之間取舍。文獻[8，15，18]將信號用傅里葉系數表示，對于帶限信號，可以只保留有限個非零系數來近似信號。因此，在這里將截獲信號用傅里葉系數表示，同時考慮量化和同步等誤差。

(6)

可以將式(6)寫成如下向量形式：

(7)

式中：

Dl,k=diag(e-j2πf-Ntl,k,e-j2πf-N+1tl,k,…,e-j2πfNtl,k)

綜合相關文獻給出輻射源位置直接估計方法的步驟如下。將所有觀測站截獲得到的輻射源信號集中傳輸到同一位置，如主觀測站。對于每一個網格點pg，g=1,2,…,G，計算信號從該網格點傳播到各觀測站的傳播時間和多普勒頻率，構造如下包含位置信息的矩陣:

(8)

式中：矩陣Dl和Fl對應著網格點傳播到各觀測站的傳播時間和多普勒頻率。

對于任意給定的矩陣X，由于XHX與XXH非零特征值相同。因此，構造L×L矩陣：

(9)

將所有時間段該矩陣最大特征值相加，將各網格點對應值組成集合，集合中的最大值對應的網格點即為輻射源估計位置?？蓪懗?/p>

(10)

矩陣Qpg,k為互模糊矩陣。Qpg,k對角線上的第l個元素為觀測站l在第k時間段內截獲信號的能量。Qpg,k的第(l,m)個非對角線元素是觀測站l和m在第k個時間段內截獲信號的互模糊函數。

(11)

式中：τg和fg分別為觀測站l和m相對于網格點pg所對應的時差和多普勒頻差。

2 去中心化直接定位方法

相關直接定位方法構造互模糊矩陣Qpg,k用于特征值求解以估計輻射源位置，都需要計算出所有矩陣元素。由于互模糊矩陣是Hermite矩陣，因此對于每個觀測時間段都需要計算L(L-1)/2個互模糊函數和L個截獲信號能量。相關直接定位方法需將各觀測站截獲信號轉發到主觀測站，再針對每個網格點將所有觀測站截獲信號構造成數據矩陣以構造互模糊矩陣。

2.1 去中心化構造僅滿秩互模糊矩陣

兩步定位方法利用觀測站截獲信號測量的時差頻差進行定位解算時，一般只利用了L-1組時差頻差。而將L個觀測站兩兩配對，最多可以有LL-1/2組配對，計算LL-1/2個互模糊函數，測量得到LL-1/2組時差頻差。直接定位方法構造互模糊矩陣Qpg,k時利用了全部LL-1/2 組配對，而兩步定位方法一般只利用其中L-1組實現定位解算。直接定位法構造的互模糊矩陣Qpg,k如果也只利用L-1組配對，計算L-1個互模糊函數構造僅滿秩的互模糊矩陣，不構造所有矩陣元素，此時互模糊矩陣剛好滿足滿秩條件，可求解特征值以估計輻射源位置。構造僅滿秩的互模糊矩陣可降低定位系統的計算量。

根據文獻[24]的思路，將時差頻差用變換矩陣乘以到達時間和到達頻率向量進行表示。令τd=Tτ和fd=Tf分別表示定位中涉及到的時差和頻差向量。其中:τ=[τ1τ2…τL]T和f=[f1f2…fL]T分別為L個觀測站截獲的信號到達時間向量和到達頻率向量；T為列數為L的變換矩陣，T的行數可以為L-1至L(L-1)/2的任意值，T的每行只有一個1和-1。

兩步定位方法提取L-1個時差頻差需設定主觀測站，其他L-1個觀測站為副觀測站，將主副觀測站配對，各副觀測站截獲信號轉發到主觀測站計算互模糊函數，提取中心化時差頻差。不失一般性假定觀測站1為主觀測站，中心化時差向量和頻差向量分別為τdC=[τ21τ31…τL1]T和fdC=[f21f31…fL1]T。中心化變換矩陣為

(12)

根據文獻[25]，與觀測站l配對的觀測站集合用Nl表示，與觀測站l配對的觀測站數目稱為觀測站l的度(degree)，用Dl表示。Dl的數值即對應著互模糊矩陣第l行和第l列中非零元素的數目。兩步定位法的中心化配對方案主觀測站的度D1=L，副觀測站的度Dl=2。

中心化配對方案計算構造的僅滿足滿秩條件的互模糊矩陣中各非對角線非零元素分別對應不同觀測站配對的互模糊函數，對角線元素為各觀測站截獲信號的能量。具體形式為

(13)

式中：All為觀測站l截獲信號的能量；Clm為觀測站l和m對應網格點pg的互模糊函數數值。

中心化時差頻差的互模糊函數計算和提取需將數據傳輸和計算都集中于單一觀測站或者融合中心，限制了定位的數據傳輸和計算效率以及觀測站網絡擴展。由于直接定位方法無需滿足中心化限制，可采用去中心化的觀測站配對方案，不設定主觀測站，各觀測站只向鄰近站轉發截獲信號，在各觀測站并行計算L-1個互模糊函數，提高直接定位的數據傳輸和計算效率，便于觀測站網絡擴展。

由于觀測站間的基線長度影響定位性能，為了利用有限的L-1個互模糊函數盡可能地提高定位精度，可以使觀測站配對間總距離盡可能大。

(14)

目標函數式(14)得到的最優觀測站標號排列Sn使得觀測站配對間總距離最大。將Sn中觀測站標號排列依次序替換為1到L。采用去中心化配對觀測站方案，將替換標號后的觀測站l和l+1 進行配對，并行地將觀測站l+1截獲信號只轉發至鄰近觀測站l。此時各觀測站的度相同，都是Dl=3。采用分布式計算，將原來互模糊函數計算集中于單一觀測站或融合中心的方式改為分散到各觀測站并行計算，縮減計算所有互模糊函數所需時間。

去中心化觀測站配對方案對應去中心化時差和頻差向量分別為τdD=[τ21τ32…τL(L-1)]T和fdD=[f21f32…fL(L-1)]T，去中心化變換矩陣為[26]

(15)

采用去中心化配對方案計算互模糊函數構造的互模糊矩陣的具體形式為

(16)

式中：All與Clm的意義與式(13)相同。

綜上，給出去中心化構造僅滿足滿秩條件的互模糊矩陣，實現直接定位的步驟如下：

步驟1根據觀測站位置的拓撲結構，按照目標函數式(14)得到各觀測站的最優標號集合。

步驟2在各觀測站計算第k時間段內截獲信號的能量El,k。

步驟3在第k時間段，任意觀測站l+1將截獲信號轉發至觀測站l，其中l=1,2，…,L-1。

步驟4在觀測站l計算CAFl(l+1),k，并將相應的時差頻差的互模糊函數數值與對應的網格點pg關聯，得到[CAFl(l+1),k]pg。

步驟5對于每一個網格點pg，構造形如式(16)的互模糊矩陣Qpg,k，其中第l個對角線元素為El,k，第l行l+1列與第l+1行l列非對角元素分別為[CAFl(l+1),k]pg與其共軛。

步驟6計算該時間段內互模糊矩陣Qpg,k的最大特征值λmax(Qpg,k)。

2.2 去中心化補全互模糊矩陣

根據2.1節構造的互模糊矩陣估計輻射源位置，雖然極大地提高了數據傳輸和計算的效率，但因為沒有充分利用所有的觀測站配對以計算互模糊矩陣的所有元素，因此無法達到原始直接定位方法的精度性能，在低信噪比下得到精度較高的定位結果。文獻[2]在提出的HC直接定位方法中給出了主副站間的互模糊函數關系公式。本節將該關系公式推廣，采用2.1節利用信號計算L-1個去中心化互模糊函數后，再去中心化的在各觀測站并行計算余下的互模糊函數，補全互模糊矩陣所有元素。

由于觀測站間需多次數據傳輸，為了降低觀測站間的數據傳輸時延，需要觀測站配對間總距離盡可能小。

因此，根據目標函數式(14)中各符號的定義，以及各觀測站位置構造如下目標函數：

(17)

目標函數式(17)得到的最優觀測站標號排列Sr使得觀測站配對間總距離最小。將Sn中的觀測站標號排列依次序替換為1到L。

因為變換矩陣T秩為L-1，而τmn=τln-τlm和fmn=fln-flm，因此在所有L(L-1)/2個時差頻差中，獲得其中包含所有觀測站數據的L-1個時差頻差，其他時差頻差都可以由L-1個時差頻差組成的向量經過初等變換計算得到[24]。

由于各觀測站截獲的是相同輻射源同一時刻的信號，在所有L(L-1)/2個互模糊函數中，計算出其中L-1個包含所有觀測站數據的互模糊函數，其他互模糊函數都可以由L-1個互模糊函數計算得到。下面將文獻[2]的關系公式推廣到任意觀測站配對的情況。第k時間段觀測站l截獲的不含噪聲信號為

(18)

式(18)即式(1)不含噪聲的形式。相同時間段內觀測站m截獲信號與式(18)的關系為

(19)

因此，兩者的互模糊函數可表示為

CAFlm,k(τ,f)=

G*ej2πflm,kτAAFl,k(τ-τlm,k,f-flm,k)

(20)

式中：

(21)

由式(21)可得

(22)

式中：AAFl,k(0,0)為觀測站l截獲信號的能量El,k。

任意觀測站m與n的互模糊函數可用任意觀測站l與n的互模糊函數表示，表達式為

CAFmn,k(τ,f)=

Gej2πfτlm,kCAFln,k(τ+τlm,k,f+flm,k)

(23)

將式(22)代入式(23)可得互模糊函數之間關系公式。

CAFmn,k(τ,f)=

(24)

利用式(24)任意互模糊函數間的關系公式，可采用歸約方式在各觀測站并行計算余下互模糊函數，補全互模糊矩陣，降低時間復雜度。根據文獻[27]，歸約(Reduction)是將位于不同觀測站的多個結果綜合起來，通過某些操作或計算產生新的結果，存儲在指定的觀測站。采用歸約的方式在各觀測站并行計算互模糊函數，共需計算L-2次循環。循環之前將觀測站l+1計算的信號能量El+1,k傳遞到觀測站l。在第i次循環中，觀測站l+1將上次循環計算得到CAF(l+1)(l+i+1),k傳遞到觀測站l，根據互模糊函數之間關系公式(24)，在觀測站l處與CAFl(l+1),k和El+1,k計算得到CAFl(l+i+1),k，再執行下次循環。第i次循環中l=1,2，…,L-i-1。

綜上，基于L-1個去中心化互模糊函數，給出利用互模糊函數之間關系，采用歸約的方式，補全互模糊矩陣，實現直接定位的步驟如下：

步驟1根據觀測站位置的拓撲結構，按照目標函數式(17)得到各觀測站的最優標號集合。

步驟2在各觀測站l計算第k時間段內截獲信號的能量El,k，并將其轉發至觀測站l-1。

步驟3和步驟4與2.1節中對應的步驟3和步驟4相同。

步驟6依次執行如下循環操作，在第i=1,2，…,L-2 次循環中，觀測站l=1,2，…,L-i-1并行執行：

(1) 接收從觀測站l+1傳遞的上次循環的計算結果CAF(l+1)(l+i+1),k。

(3) 將相應的時差頻差的CAFl(l+i+1),k數值與對應的網格點pg關聯，得到[CAFl(l+i+1),k]pg。

步驟7對于每一個網格點pg，構造完整的互模糊矩陣Qpg,k，其中第l個對角線元素為El,k，第l行m列與第m行l列非對角元素分別為[CAFlm,k]pg與其共軛。

步驟8和步驟9與2.1節中步驟6和步驟7相同。

3 定位性能分析

下面分析比較定位方法的誤差克拉美羅下界與計算復雜度。

3.1 定位誤差克拉美羅下界

克拉美羅下界可以通過Fisher信息矩陣(Fisher Information Matrix, FIM)的逆求出[28]。為了分析的簡便，對二維空間情況進行推導分析，三維空間的情況可以經過拓展得到。

文獻[17]給出的觀測站l截獲信號經過變換和量化等之后傳遞到觀測站m，任意時間段k內采用互模糊函數方法估計時差頻差的FIM如下：

Jlm=2Re{GH[Σ+W]-1G}

(25)

式中：Σ為零均值加性高斯白噪聲協方差；W為量化誤差協方差；

(26)

為了便于分析，假定各觀測站配對用互模糊函數方法估計的FIM相同，且滿足：

(27)

式中：結合文獻[6，22]有

(28)

(29)

根據文獻[17，26]，相關直接定位方法所有觀測站兩兩配對計算互模糊函數構造互模糊矩陣的所有元素，輻射源位置估計的FIM為

(30)

式中：

(31)

當l≠m時，Tlm為變換矩陣T中某一行，該行的第l列為1，第m列為-1，其他都為零；當l=m時，行向量Tll中所有元素都為零。

將式(27)代入式(30)，整理得到FIM各元素為

(32)

同理，兩步定位方法中采用中心化觀測站配對方案，輻射源位置估計的FIM為

(33)

整理得到兩步定位方法的FIM各元素為

(34)

2.1節計算L-1個互模糊函數去中心化的構造僅滿秩互模糊矩陣的輻射源位置估計FIM如下：

(35)

整理得到該FIM各元素為

(36)

(37)

(38)

(39)

比較式(37)～式(39) 3個公式，相關直接定位方法各觀測站配對都充分利用，FIM對角線上元素數值最大，因此位置估計的CRLB最小，精度性能最穩定。中心化配對方案的FIM非零元素集中在第1行和第1列，每個觀測站配對都利用到主觀測站與輻射源在空間中的相對位置、速度和截獲信號等信息，位置估計CRLBC受主觀測站影響較大，精度性能不穩定。去中心化的FIM非零元素均勻分布，均勻利用各觀測站與輻射源在空間中相對位置、速度和截獲信號，隨著輻射源在空間中不同位置的變化，位置估計CRLBD具有較穩定的起伏，甚至可達到更高的位置估計精度。但與相關直接定位方法CRLB相比，精度性能稍差。

3.2 計算復雜度

中心化方法將所有計算集中于單一觀測站或融合中心，而去中心化方法將計算分散到各觀測站。觀測站的數據傳輸量和計算量與觀測站的度Dl有關，原始直接定位法和兩步定位方法中主觀測站的度D1=L，而本文所提兩種方法各觀測站度Dl=3。

原始直接定位法在主觀測站或融合中心需要接收所有L-1個觀測站的截獲信號。對于傅里葉系數長度為2N的信號，構造互模糊矩陣所需的L個觀測站的截獲信號能量和L(L-1)/2個互模糊函數數值的計算量為O2N，互模糊矩陣特征值分解的復雜度為OL3。因此對于格點數目為G的空域，總計算量為O([GL2-(G-2)L]·N+GL3)，且所有計算量都集中在主觀測站或融合中心。

本文第1種方法各觀測站只需接收一個觀測站的截獲信號，另外將G維的互模糊函數值傳輸到融合中心。每個觀測站只需計算自身截獲信號能量和一個互模糊函數，每個觀測站的總計算量為O2(G+1)N。L-1個觀測站并加上特征值分解的總計算量為O2(L-1)(G+1)N+GL3。

第2種方法增加的總計算量是在第1種方法基礎之上再利用式(24)互模糊函數之間關系公式計算余下(L-1)(L-2)/2個互模糊函數，而此時每計算一個互模糊函數的計算量為OG。因此與第1種方法相比增加計算量為O(0.5(L2-3L+2)G)，第2種方法總計算量為O(2(L-1)·(G+1)N+0.5(2L3+L2-3L+2)G)。

4 仿真實驗

采用均方根誤差(Root Mean Squares Error，RMSE)度量輻射源位置的定位估計精度性能，其定義為

(40)

下面對各算法進行仿真實驗分析。仿真中假定信號樣式為線性調頻信號，載頻3 GHz，調頻帶寬10 MHz，脈沖寬度1 ms，脈沖間隔10 ms。各觀測站的接收機噪聲帶寬都為20 MHz，信噪比相等。采用文獻[31]的兩步定位方法(Two-step)和文獻[7]的原始直接定位方法(Original DPD)與本文所提兩種方法進行仿真比較。本文所提兩種方法依次分別為：第1種方法(Proposed method 1)采用去中心化配對方案計算L-1個互模糊函數構造的僅滿秩互模糊矩陣的直接定位方法，第2種方法(Proposed method 2)利用任意互模糊函數間關系公式采用歸約方式，去中心化的在各觀測站并行計算余下互模糊函數以補全互模糊矩陣的直接定位方法。觀測站與輻射源的位置與速度如表1所示。

在該觀測站位置拓撲結構下，如果各觀測站之間都可以單跳傳輸。采用去中心化配對方案，遍歷各觀測站配對集合。本文所提第1種方法中的目標函數式(14)的最優配對方案為{(1,4)，(4,7)，(7,9)，(9,5)，(5,3)，(3,2)，(2,6)，(6,8)}，總傳輸距離為27 713 m，因此對本文所提第1種定位方法進行仿真時采用該最優配對方案，對觀測站標號集合{1,4,7,9,5,3,2,6,8}重新依次編號。本文所提第2種方法中的目標函數式(17)的最優配對方案為{(2,8)，(8,3)，(3,9)，(9,4)，(4,6),(6,5),(5,1),(1,7)}，初次循環時總傳輸距離為10 390 m，因此對本文所提第2種定位方法進行仿真時采用該最優配對方案，對觀測站標號集合{2,8,3,9,4,6,5,1,7}重新依次編號。

表1 觀測站和輻射源的位置與速度Table 1 Position and velocity of sensors and emitters

兩步定位方法采用的中心化觀測站配對方案，對應的CRLBC根據式(34)的Fisher信息矩陣計算。本文所提第1種方法采用的去中心化觀測站配對方案，對應的CRLBD根據式(36)的Fisher信息矩陣計算。本文所提第2種方法和原始直接定位方法采用的所有觀測站兩兩配對方案，對應的CRLB根據式(32)的Fisher信息矩陣計算。對單個輻射源進行仿真，依次對近距離和遠距離輻射源進行定位估計，各方法的RMSE和各觀測站配對方案的CRLB隨著信噪比的變化曲線如圖1所示。

圖1 各方法下不同輻射源定位隨信噪比(SNR)變化的RMSE和CRLB曲線 Fig.1 Curves of RMSE of different methods and CRLB vs signal-to-noise ratio (SNR) for different emitters

由圖1還可以看出，在相同的信噪比情況下，采用去中心化配對方案的本文所提第1種方法由于均勻利用各觀測站與輻射源在空間中相對位置，精度性能優于采用中心化配對方案的兩步定位方法。本文所提第2種方法和原始直接定位方法由于采用所有觀測站兩兩配對方案，利用了L個觀測站所有LL-1/2組兩兩配對，精度明顯優于只利用L個觀測站中L-1組兩兩配對的兩步定位方法和本文所提第1種方法；由于估計的去中心化時差頻差存在誤差，影響采用式(24)計算的余下互模糊函數的精度，使得最終的定位精度沒有達到原始直接定位方法的精度，但計算量和數據傳輸量與原始直接定位方法相比得到了極大降低。

下面仿真不同觀測站配對方案對單個輻射源的定位誤差幾何分布。由于輻射源目標在空間中的不同位置引起的定位誤差是不同的，而定位誤差幾何分布描述了不同空間位置定位誤差的分布情況。在仿真中，假定各觀測站接收機的信噪比都為0 dB。各觀測站配對方案的定位誤差幾何分布分別如圖2所示。

由圖2可以看出，由于各觀測站位置都相同，因此定位誤差在空間位置中的分布趨勢相似。與前面RMSE仿真分析結論相同，去中心化觀測站配對方案的定位精度高于中心化觀測站配對方案的定位精度，所有觀測站兩兩配對的方案在空間中的定位精度最高。

圖2 定位誤差幾何分布 Fig.2 Geometric distribution of location errors

仿真實驗采用大時寬帶寬的線性調頻脈沖信號，時頻耦合效應會使得互模糊函數有較高的副峰。當信噪比較低時，時頻耦合效應引起的較高副峰與噪聲形成的峰相加會出現超過主峰的情況，使得真實峰值模糊，無法得到正確的定位估計結果。在本文的仿真參數下，本文所提方法在信噪比低于-30 dB時會出現定位結果錯誤的情況，而原始直接定位方法在信噪比低于-35 dB時才會出現定位結果錯誤的情況?？赏ㄟ^加窗抑制副峰，改善時頻耦合效應對定位結果的影響。但加窗抑制副峰是以增大主峰寬度為代價的，會降低定位精度性能。

另外，同頻多輻射源問題一直以來都是實現無源定位的難題，互模糊函數會出現偽峰，無法得到正確的定位估計結果，尤其是當同頻多輻射源相距較近或信號功率相差較大時。對同頻輻射源可利用信號其他域特征的差異將同頻多目標的信號進行關聯，具體信號特征包括到達方向、到達時間段、調制樣式和極化等。

5 結 論

針對原始利用時差頻差的直接定位方法存在數據傳輸量和計算量大的瓶頸，本文提出了兩種去中心化直接定位方法，不同程度地簡化了構造互模糊矩陣的方法，從而降低總的計算量，并將數據傳輸和計算量分散到各觀測站之間，而不是集中于主觀測站。

1) 第1種方法不需計算全部互模糊矩陣元素，L個觀測站只計算L-1個互模糊函數，構造僅滿足滿秩條件的互模糊矩陣，采用去中心化配對方案，只將各觀測站截獲信號在站間進行一次傳輸，將數據傳輸和計算分散到各觀測站間并行執行，提高了直接定位數據傳輸和計算的效率。

2) 第2種方法在利用信號計算L-1個去中心化互模糊函數之后，根據推導的任意互模糊函數間關系公式，采用歸約方式去中心化的在各觀測站并行計算余下的互模糊函數，補全互模糊矩陣所有元素，降低了直接定位數據傳輸和計算量。

3) 性能分析和仿真實驗表明本文兩種去中心化直接定位方法精度性能優于兩步定位方法；在低信噪比時兩種方法可得到比較理想的精度性能；在高信噪比時第2種方法與原始直接定位方法的精度性能相當。

[1] 趙國慶. 雷達對抗原理[M]. 第2版. 西安: 西安電子科技大學出版社, 2012: 9-11.

ZHAO G Q. Principle of radar countermeasure[M]. 2nd ed. Xi’an: Xidian University Press, 2012: 9-11 (in Chinese).

[2] POURHOMAYOUN M, FOWLER M L. Distributed computation for direct position determination emitter location[J]. IEEE Transactions on Aerospace and Electronic Systems, 2014, 50(4): 2878-2889.

[3] WAX M, KAILATH T. Decentralized processing in sensor arrays[J]. IEEE Transactions on Acoustics, Speech, and Signal Processing, 1985, 33(4): 1123-1129.

[4] WEISS A J. Direct position determination of narrowband radio frequency transmitters[J]. IEEE Signal Processing Letters, 2004, 11(5): 513-516.

[5] WEISS A J. Direct position determination in the presence of model errors-known waveforms[J]. Elsevier Digital Signal Processing, 2006, 16(1): 52-83.

[6] AMAR A, WEISS A J. Localization of radio emitters based on Doopler frequency shifts[J]. IEEE Transactions on Signal Processing, 2008, 56(11): 5500-5508.

[7] WEISS A J. Direct geolocation of wideband emitters based on delay and Doopler[J]. IEEE Transactions on Signal Processing, 2011, 59(6): 2513-2521.

[8] OISPUU M, NICKEL U. Direct detection and position determination of multiple sources with intermittent emission[J]. Signal Processing, 2010, 90: 3056-3064.

[9] 張敏, 郭福成, 周一宇. 基于單個長基線干涉儀的運動單站直接定位[J]. 航空學報, 2013, 34(2): 378-386.

ZHANG M, GUO F C, ZHOU Y Y. A single moving observer direct position determination method using a long baseline interferometer[J]. Acta Aeronautica et Astronautica Sinica, 2013, 34(2): 378-386 (in Chinese).

[10] 張敏, 郭福成, 周一宇, 等. 運動單站干涉儀相位差直接定位方法[J]. 航空學報, 2013, 34(9): 2185-2193.

ZHANG M, GUO F C, ZHOU Y Y, et al. A single moving observer direct position determination method using interferometer phase difference[J]. Acta Aeronautica et Astronautica Sinica, 2013, 34(9): 2185-2193 (in Chinese).

[11] 王鼎, 張剛, 沈彩霞, 等. 一種針對恒模信號的運動單站直接定位算法[J]. 航空學報, 2016, 37(5): 1622-1633.

WANG D, ZHANG G, SHEN C X, et al. A direct position determination algorithm for constant modulus signals[J]. Acta Aeronautica et Astronautica Sinica, 2016, 37(5): 1622-1633 (in Chinese).

[12] BOSSE J, FERREOL A, GERMOND C, et al. Passive geolocalization of radio transmitters: Algorithm and performance in narrow band context[J]. Signal Processing, 2012, 92: 841-852.

[13] 黃志英, 吳江, 唐濤, 等. 非一致噪聲下的多陣列直接定位算法[J]. 西安交通大學學報, 2015, 49(10): 136-142.

HUANG Z Y, WU J, TANG T, et al. A direct position determination algorithm based on multi arrays in the presence of unknown nonuniform noise[J]. Journal of Xi’an Jiaotong University, 2015, 49(10): 136-142 (in Chinese).

[14] VANKAYALAPATI N, KAY S M. Asymptotically optimal localization of an emitter of low probability of intercept signals using distributed sensors[J]. IEEE Transactions on Aerospace and Electronic Systems, 2012, 48(1): 737-748.

[15] VANKAYALAPATI N, KAY S M, DING Q. TDOA based direct positioning maximum likelihood estimator and the Cramer-Rao bound[J]. IEEE Transactions on Aerospace and Electronic Systems, 2014, 50(3): 1616-1635.

[16] FOWLER M L, CHEN M, BINGHAMTON S. Fisher-information-based data compression for estimation using two sensors[J]. IEEE Transactions on Aerospace and Electronic Systems, 2005, 41(3): 1131-1137.

[17] CHEN M, FOWLER M L. Data compression for multiple parameter estimation with application to emitter location systems[J]. IEEE Transactions on Aerospace and Electronic Systems, 2010, 46(1): 308-322.

[18] 王云龍, 吳瑛. 聯合時延與多普勒頻率的直接定位改進算法[J]. 西安交通大學學報, 2015, 49(4): 123-129.

WANG Y L, WU Y. An improved direct position determination algorithm with combined time delay and Doppler[J]. Journal of Xi’an Jiaotong University, 2015, 49(4): 123-129 (in Chinese).

[19] XIA W, LIU W. Distributed adaptive direct position determination of emitters in sensor networks[J]. Signal Processing, 2016, 123: 100-111.

[20] XIA W, LIU W, ZHU L F. Distributed adaptive direct position determination based on diffusion framework[J]. Journal of Systems Engineering and Electronics, 2016, 27(1): 28-38.

[21] KNAPP C H, CARTER G C. The generalized correlation method for estimation of time delay[J]. IEEE Transactions on Acoustics, Speech, and Signal Processing, 1976, 24(4): 320-327.

[22] STEIN S. Algorithms for ambiguity function processing[J]. IEEE Transactions on Acoustics, Speech, and Signal Processing, 1981, 29(3): 588-599.

[23] 朱穎童, 董春曦, 劉松楊, 等. 存在觀測站位置誤差的轉發式時差無源定位[J]. 航空學報, 2016, 37(2): 706-716.

ZHU Y T, DONG C X, LIU S Y, et al. Passive localization using retransmitted tdoa measurements in the presence of sensor position errors[J]. Acta Aeronautica et Astronautica Sinica, 2016, 37(2): 706-716 (in Chinese).

[24] MENG W, XIE L H, XIAO W D. Decentralized TDOA sensor pairing in multihop wireless sensor networks[J]. IEEE Signal Processing Letters, 2013, 20(2): 181-184.

[25] CATTIVELLI F S, SAYED A H. Diffusion LMS strategies for distributed estimation[J]. IEEE Transactions on Signal Processing, 2010, 58(3): 1035-1048.

[26] TORRIERI D J. Statistical theory of passive location systems[J]. IEEE Transactions on Aerospace and Electronic Systems, 1984, 20(2): 183-198.

[27] 張林波. 并行計算導論[M]. 北京: 清華大學出版社, 2006: 52-53.

ZHANG L B. Introduction to parallel computing[M]. Beijing: Tsinghua University Press, 2006: 52-53 (in Chinese).

[28] KAY S M. Fundamentals of statistical signal processing: Estimation theory[M]. Upper Saddle River, NJ: Prentice Hall, 1993: 39-40.

[29] 郭福成, 樊昀, 周一宇, 等. 空間電子偵察定位原理[M]. 北京: 國防工業出版社, 2012: 155-158.

GUO F C, FAN Y, ZHOU Y Y, et al. Localization principles in space electronic reconnaissance[M]. Beijing: National Defense Press, 2012: 155-158 (in Chinese).

[30] AMAR A, LEUS G, FRIEDLANDER B. Emitter localization given time delay and frequency shift measurements[J]. IEEE Transactions on Aerospace and Electronic Systems, 2012, 48(2): 1826-1837.

[31] HO K C, XU W W. An accurate algebraic solution for moving source location using TDOA and FDOA measurements[J]. IEEE Transactions on Signal Processing, 2004, 52(9): 2453-2463.

(責任編輯： 蘇磊)

URL：www.cnki.net/kcms/detail/11.1929.V.20170109.1550.002.html

DecentralizeddirectpositiondeterminationmethodbasedonTDOAandFDOA

ZHUYingtong*,DONGChunxi,DONGYangyang,XUJin,ZHAOGuoqing

KeyLaboratoryofElectronicInformationCountermeasureandSimulationTechnology,MinistryofEducation,XidianUniversity,Xi’an710071,China

Toovercomethebottleneckofdatatransmissioncapacityandcomputationoftheoriginaldirectpositiondetermination(DPD)methodbasedonTDOAandFDOA,twodecentralizedDPDmethodsareproposed.Thefirstproposedmethodusesthedecentralizedsensorcouplingscheme,inwhichthesignalinterceptedbyeachsensoristransmittedonlyonceanddatatransmissionandcomputationoperationsaredispersedtoeachsensortocomputecrossambiguityfunction(CAF)inparallel,soastoconstructcrossambiguitymatrix(CAM)withmerelyfullrank.BasedonthederivationofaformulafortherelationbetweenarbitraryCAFs,thesecondmethodappliesareductionoperationtocomputetheremainingCAFsateachsensorinparallel,andcomplementsallelementsofCAM.TwodecentralizedDPDmethodsreducetheamountofdatatransmissionandimprovetheefficiencyofcomputation.Performanceanalysisandsimulationresultsshowthattheaccuracyoftheproposedmethodsaresuperiortothatofthetwo-stepmethods.AtlowSNR,bothofthetwoproposedmethodscanachievedesiredaccuracy,andathighSNR,thesecondproposedmethodcanobtaintheaccuracysimilartothatoftheoriginalDPDmethod.

directpositiondetermination;passivelocalization;crossambiguityfunction;decentralized;sensorcoupling;Cramer-Raolowerbound

2016-08-30；Revised2016-11-28；Accepted2016-12-22；Publishedonline2017-01-091550

s：NationalBasicResearchProgramofChina(61**81);NationalHigh-techResearchandDevelopmentProgramofChina(2014AA80**086H);theFundamentalResearchFundsfortheCentralUniversities(JB140203)

.E-mailzhuyt_xd@163.com

2016-08-30；退修日期2016-11-28；錄用日期2016-12-22； < class="emphasis_bold">網絡出版時間

時間：2017-01-091550

www.cnki.net/kcms/detail/11.1929.V.20170109.1550.002.html

國家“973”計劃 (61**81)； 國家“863”計劃 (2014AA80**086H)； 中央高?；究蒲袠I務費專項資金 (JB140203)

.E-mailzhuyt_xd@163.com

朱穎童， 董春曦， 董陽陽， 等． 去中心化時差頻差直接定位方法J．航空學報,2017,38(5):320727.ZHUYT,DONGCX,DONGYY,etal.DecentralizeddirectpositiondeterminationmethodbasedonTDOAandFDOAJ.ActaAeronauticaetAstronauticaSinica,2017,38(5):320727.

http://hkxb.buaa.edu.cnhkxb@buaa.edu.cn

10.7527/S1000-6893.2016.320727

V247.5； TN971

A

1000-6893(2017)05-320727-13