基于稀疏互質L型陣列的二維測向算法

邵 華 蘇衛民 顧 紅 王 燦

（南京理工大學電子工程與光電技術學院，江蘇 南京 210094）

引 言

在雷達、聲納和通信等領域，對輻射源進行定位［1－3］，常需要先確定入射信號的二維波達角（2DDOA，包括方位角和俯仰角）。針對2D－DOA估計問題，陣列信號處理中已有許多行之有效的方法被提出，例如二維多重信號分類法［4］（2D－MUSIC），但大多數測向方法只適合于等間距面陣。與等間距面陣相比，稀疏面陣能在不損失陣列孔徑的前提下減小所需的陣元數，因此，研究稀疏面陣結構及其相應的波達角（DOA）估計算法具有極大的意義。

其中，文獻［5］提出二維旋轉不變（2D－ESPRIT）算法，它通過子陣間x軸和y軸的旋轉不變關系估計方位角和俯仰角，缺點是需要3個結構相同的子陣，硬件成本高；基于此缺點，文獻［6］、［7］利用四階累積量［8－10］的陣列擴展特性，分別沿x軸和y軸構造與原陣列結構相同的虛擬陣列，然后再由2D－ESPRIT算法進行2D－DOA估計，即二維虛擬旋轉不變（2D－VESPA）算法，然而其要求參考陣元間距（原陣列參考陣元和虛擬陣列參考陣元之間的間距）不大于信號半波長0.5λ；為突破陣元間距不大于0.5λ的限制條件，文獻［11］、［12］提出嵌套陣列，文獻［13］、［14］提出互質陣列，利用子陣間的嵌套或互質關系形成陣元間距不大于0.5λ的均勻矩形虛擬陣列，而允許物理陣元間距大于0.5λ.雖然這兩種算法是利用虛擬陣元消除了測向模糊，但是其實質是僅采用一步就完成了解模糊。這導致正確解模糊時，系統允許兩個陣元接收信號間相位差的測量誤差（簡稱系統容差）較小。

針對以上問題，本文基于模轉換［15］逐步解模糊的思想，首先構造了一類由兩個具有相同互質結構的線陣組成的L型陣列，其中線陣1沿x軸放置，線陣2沿y軸放置；其次，利用一維虛擬旋轉不變（1DVESPA）算法估計各基線所對應的模糊相位差；然后，結合基線間的互質關系逐步解模糊，分別得到x軸和y軸中最長基線對應的無模糊方向余弦；最后，采用x軸和y軸特征向量（估計不同坐標軸方向余弦時分別進行特征分解）間的對應關系，實現不同坐標軸方向余弦的匹配，最終得到方位角和俯仰角的估計值。與2D－VEPSA算法相比，本算法最大的優勢在于利用基線間的互質關系突破0.5λ的限制條件，提高測向精度的同時保證較大的系統容差。

1 信號模型

如圖1所示，L型陣列由x軸子陣（M＋1個全向陣元）和y軸子陣（M＋1個全向陣元）組成，兩個子陣在坐標原點處共用陣元0，并且兩子陣內相鄰陣元間的基線長度分別為Dx，1，…，Dx，M和Dy，1，…，Dy，M，其中Dx，m＝Dy，m＝Dm，m＝1，…，M，基線D1最長，各基線之間滿足如下比例關系

式中，Pm和Qm為互質的正整數，并且Q2＜…＜QM.

圖1 L型陣列結構示意圖

假設K個入射方向為（θ1，φ1），…，（θK，φK）的非高斯獨立窄帶遠場信號源入射到上述陣元總數為2M＋1的L型陣列天線上，其中0≤φk＜2π和0≤θk＜π分別表示第k個入射信號的方位角和俯仰角。那么，L型陣列在t時刻接收的信號矢量為

式中：yx，m（t）和yy，m（t）分別表示在x軸和y軸上第m個陣元的接收信號；s（t）＝ ［s1（t），…，sK（t）］T為信號向量；w（t）＝ ［w0（t），…，w2M（t）］T表示陣列接收噪聲，噪聲為平穩、時間和空間都互不相關的高斯白噪聲，且與信號相互獨立；A＝ ［a1，…，aK］是（2M＋1）×K維的陣列導向矩陣；ak＝］T，k＝1，…，K，其中

式中：qx，k，m＝exp（－jψx，k，m）和qy，k，m＝exp（－jψy，k，m），m＝1，…，M，其中ψx，k，m＝／λ和＝2πvkDy，m／λ，uk＝sinθkcosφk和vk＝sinθksinφk分別表示x軸和y軸方向余弦。當基線Dx，m＝Dy，m＞λ／2時，測量相位差φx，k，m∈ ［0，2π）和φy，k，m∈ ［0，2π）具有周期性模糊，即

式中，Kx，k，m∈Z和Ky，k，m∈Z分別表示φx，k，m和φy，k，m的模糊數。這導致x軸和y軸測量方向余弦也具有模糊性。利用L型陣列輸出的N個快拍數據進行無模糊2D－DOA估計。另外，為表述簡潔，下文均省略時間t.

2 算法描述

2.1 四階累積量矩陣的構造

利用圖1中L型陣列的特殊結構，并結合四階累積量的定義［7－10］及空間接收信號的數學模型，構造如下四階累積量矩陣

式中，m＝1，…，M，且

式中，γ4，sk表示第k個信號的四階累積量。由四階累積量的陣列擴展特性可知，式（8）表示的物理含義為：原陣列沿x軸平移xm距離后產生的虛擬陣列，如圖2（a）所示，其中陣元m為虛擬陣列的參考陣元，稱這類虛擬陣列為x軸虛擬陣列組；式（9）則表示原陣列沿y軸平移ym距離后產生的虛擬陣列，如圖2（b）所示，其中陣元M＋m為虛擬陣列的參考陣元，類似的稱其為y軸虛擬陣列組。為了估計所有相鄰陣元間的相位差，選擇陣元0作為原陣列的參考陣元，其他2M個陣元分別為2M個虛擬陣列的參考陣元，那么可以定義如下矩陣

對式（13）進行奇異值分解（SVD）后，可得信號子空間Es＝ ［，…，，…，］T.

2.2 估計x軸基線對應的測量相位差

抽取A的相應行可以構造原陣列和x軸虛擬陣列組對應的導向矩陣Ax＝ ［AT，ATΛx，1，…，ATΛx，M］T.同樣由Es可得它們對應的信號子空間Es，x＝ ［，，…，］T.因為Es，x和Ax都表示信號子空間，所以它們之間可以相互轉換，令可逆矩陣T為轉換矩陣，則Es

為了估計基線對應的測量相位差，令Fx，m＝（Es，x，m－1）?Es，x，m，m＝1，…，M，其中Es，x，0＝Es，0，并結合式（16）得

其中：（·）?表示矩陣偽逆；Λx，0是元素為0的K維方陣。對Fs，x，m進行特征分解后得Λx，m－Λx，m－1，求其對角元素的相位可得測量相位差的估計值＝，…］和特征向量Tx，m.

2.3 x軸方向余弦的估計

由于基線越長，測向精度越高。因此，為得到高精度x軸方向余弦估計，必須先通過模轉換解模糊算法估計最長基線Dx，1對應的無模糊相位差。然而，由模轉換解模糊算法的原理可知，它需要結合相同信號、不同基線所對應的測量相位差才能正確解模糊。又由x軸基線對應的測量相位差的求解過程可知，和是分別進行特征分解后得到的，其中m≠n，m∈1［］M，n∈［1M］，故和信號源三者之間的順序具有任意性，所以對它們的配對是正確解模糊的前提。

根據矩陣的特征分解可知，特征值與特征向量是一一對應的，即和的配對問題轉化成Tx，m和Tx，n的配對問題。假設表示Tx，m中第k1列向量表示Tx，n中第k2列向量，當它們對應 相 同 信 號 源 時，＝ 1，否 則?1.基于以上分析，對Tx，m和Tx，n做如下操作其中，矩陣G中1的位置反映和間的配對關系。不失一般性，假設以的順序為參考，利用上述方法分別對｛，m＝1，…，M｝進行配對后得

結合基線Dx，m和配對后的測量相位差估計值，m＝1，…，M，k＝1，…，K，利用模轉換解模糊算法解最長基線D1，x對應測量相位差估計值φx，k，1的模糊，得到Dx，1對應的無模糊測量相位差，再經簡單的計算后得高精度、無模糊x軸方向余弦估計值經過逐步解模糊得出，與單步解模糊相比，它能保證較大的系統容差較［15］。

2.4 y軸方向余弦的估計

與x軸方向余弦的估計類似，y軸方向余弦的估計分為以下幾步：

1）由Es和A分別構造原陣列和y軸虛擬陣列組對應的信號子空間Es，y＝ ［，，…，］T和導向矩陣Ay＝ ［AT，ATΛy，1，…，ATΛy，M］T；

2）為了估計基線對應的測量相位差，m＝1，…，M，令

其中，Es，y，0＝Es，0，Λy，0是元素為0的K維方陣；

3）對Fs，y，m進行特征分解后得基線Dy，m對應測量相位差的估計＝［，…，］和特征向量Ty，m，m＝1，…，M；

4）以為基準，利用G＝，m，m＝2，…，M中1的位置對和進行配對，最后得到配對后y軸各基線所對應的測量相位差

5）對于第k個信號，利用模轉換算法解模糊相位差的模糊數，得到無模糊測量相位差再經簡單的計算后得高精度、無模糊y軸方向余弦估計

2.5 方位角和俯仰角的估計

由于和間并不一定對應同一信號源，因此，先對它們進行配對。與前面的配對類似，采用G＝中1的位置實現和的匹配，記匹配后的x軸和y軸方向余弦分別為和.根據方向余弦的定義，可得第k個入射信號方位角和俯仰角的估計值

3 計算機仿真試驗

通過計算機仿真來比較本文算法和2D－VEPSA算法的測向性能。由于本文算法要求各陣元間距滿足一定的互質關系，而2D－VEPSA算法要求參考陣元間距不大于信號半波長，其他陣元位置任意，即它們對陣列結構的要求不同。因此，在仿真中，兩種算法分別采用陣元數相同、結構不同的天線a和天線b來比較它們的測向性能，其中天線a是陣元間距滿足互質關系的9陣元L型陣，各陣元都為全向陣，并且為了減小天線a的體積，已對輔助陣元的位置做了調整，但并不影響基線間的比例關系，如圖3所示；天線b也是9陣元L型陣，其陣元1和陣元5的位置為 （λ／2，0，0）和（0，λ／2，0），其他陣元分布與天線a的相同，它作為天線a對比天線。

圖3 仿真中本文算法采用的天線a

必須說明的是：在仿真中，本文算法首先以陣元0到陣元｛m，m＝1，…，4｝為參考陣元，形成與原陣列相距 ｛（Dm，0，0），m＝1，…，4｝的x軸虛擬陣列組，然后以陣元0到陣元｛m，m＝5，…，8｝為參考陣元，形成與原陣列相距 ｛（0，Dm，0），m＝1，…，4｝的y軸虛擬陣列組，利用虛擬陣列與原陣列間的旋轉不變關系估計兩坐標軸各基線對應的模糊相位差，經解模糊得到兩坐標軸最長基線（D1）對應的無模糊方向余弦，簡單計算后得三個信號的2D－DOA估計值；2D－VEPSA算法以陣元0、1和5為參考陣元，形成與原陣列相距 （λ／2，0，0）和（0，λ／2，0）的虛擬陣列組，利用虛擬陣列組與原陣列的旋轉不變關系估計陣元0和陣元1及陣元0和陣元5之間對應的無模糊方向余弦，從而估計三個信號的2D－DOA.另外，仿真實驗中使用均方根誤差（RMSE）對算法測向性能進行評估，其中均方根誤差定義為

式中：和分別表示第k個目標俯仰角和方位角估計值；N表示Monte－Carlo試驗次數。

仿真A：假設存在3個非高斯、獨立、窄帶遠場信號源投射到這兩個天線上，3個信源方位角分別為10°、20°和30°，俯仰角分別為15°、25°和35°，信源間相互獨立，且與噪聲互不相關，噪聲為平穩、時間和空間都不相關的高斯白噪聲。對于測量噪聲，三個信號源有相同的信噪比（SNR），SNR取值范圍為－10～20dB，取值步長為5dB，每次獨立仿真采用2 000次快拍數據，獨立重復200次Monte－Carlo試驗。

如圖4給出了兩種算法的2D－DOA估計的RMSE隨SNR的變化曲線。圖5是本文算法中正確解模糊的概率隨SNR的變化曲線。結合圖4和圖5可以看出，在正確解模糊的情況下（RSN≥0 dB），本文算法相對于2D－VEPSA算法在測向精度上具有較大改善，這主要是因為本文算法利用解模糊算法突破了2D－VEPSA算法中參考陣元間距不能大于半波長的限制，而測向精度與參考陣元的間距成正比，因此，本文算法能達到更高的測向精度。但是，當RSN＜0dB，本文算法的測向精度逐漸劣于2D－VEPSA算法，其原因是更低的信噪比超出了正確解模糊的系統容差，從而存在解模糊錯誤，導致本文算法測向性能急劇下降。

仿真B：假設兩個陣列接收信號的SNR為10 dB，當快拍數從200變化到2 000，變化步長為200時，估計兩種算法的測向性能隨快拍數的變化，其他條件與仿真A相同。

如圖6是兩種算法的2D－DOA估計的RMSE隨快拍數的變化曲線。很明顯，兩種算法的測向性能隨快拍數的增大而提高，但變化的陡峭程度不如SNR.另外，本文算法測向 RMSE要小于2DVEPSA算法，其原因也是本文算法構造稀疏互質L型陣列，利用陣元間的互質關系解模糊，擴大陣列的孔徑，提高測向精度。

圖4 兩種算法的2D－DOA估計的RMSE隨SNR的變化曲線

圖5 本文算法中正確解模糊概率隨SNR的變化

圖6 兩種算法的2D－DOA估計的RMSE隨快拍數的變化曲線

4 結 論

提出了一種基于稀疏互質L型陣列的2D－DOA估計算法。與2D－VEPSA算法相比，利用陣元間的互質關系，解決測向精度與測向模糊之間的矛盾，擴大陣列孔徑，提高測向精度，同時利用逐步解模糊，保證了較大的系統容差。

［1］DOGANCAY K. Relationship between geometric translations and TLS estimation bias in bearings－only target localization［J］.IEEE Transactions on Signal Processing，2008，56（3）：1005－1017.

［2］DOGANCAY K.Exploiting geometric translations in TLS based robot localization from landmark bearings［C］／／17th European Signal Processing Conference.Glasgow，Scotland，24－28August，2009：95－99.

［3］RAO S K，RAJESWARI K R，LINGAMURTY K S.Unscented Kalman filter with application to bearingsonly target tracking［J］.IETE Journal of Research，2009，55（2）：63－67.

［4］CHAN A Y J，LITVA J.MUSIC and maximum techniques on two－dimensional DOA estimation with uniform circular array［J］.IEE Proceedings Radar，Sonar and Navigation，1995，142（3）：105－114.

［5］KEDIA V S，CHANDNA B.A new algorithm for 2－D DOA estimation［J］.Signal Process，1997，60（3）：325－332.

［6］LIU T H，MENDEL J M.Azimuth and elevation direction finding using arbitrary array geometries［J］.IEEE Transactions on Signal Processing，1998，46（7）：2061－2065.

［7］DOGAN M C，MENDEL J M.Applications of cumulants to array processing－part I：aperture extension and array calibration［J］.IEEE Transactions on Signal Processing，1995，43（5）：1200－1216.

［8］齊子森，郭 英，王布宏，等.基于四階累積量的共形陣列波達方向估計算法［J］.電波科學學報，2011，26（4）：735－744.

QI Zisen，GUO Ying，WANG Buhong，et al.DOA estimation algorithm for conformal array based on fourthorder cumulants［J］.Chinese Journal of Radio Science，2011，26（4）：735－744.（in Chinese）

［9］沈怡平，滕升華.基于四階累積量的穩健盲波束形成算法［J］.電波科學學報，2008，23（6）：1056－1060.

SHEN Yiping，TENG Shenghua.Fourth－order cumulant based on robust blind beamforming algorithm［J］.Chinese Journal of Radio Science，2008，23（6）：1056－1060.（in Chinese）

［10］劉學斌，韋 崗，季 飛.基于四階累積量擴展孔徑的線陣設計［J］.電波科學學報，2006，21（1）：126－130.

LIU Xuebin，WEI Gang，JI Fei.Design of linear array with augmented aperture based on fourth－order cumulant［J］.Chinese Journal of Radio Science，2006，21（1）：126－130.（in Chinese）

［11］PAL P，VAIDYANATHAN P P.Nested arrays：a novel approach to array processing with enhanced degrees of freedom［J］.IEEE Transactions on Signal Processing，2010，58（8）：4167－4181.

［12］PAL P，VAIDYANATHAN P P.Two dimensional nested arrays on lattices［C］／／IEEE ICASSP.Prague，Czech Republic，May 22－27，2011：2548－2551.

［13］PAL P，VAIDYANATHAN P P.Sparse sensing with co－prime samplers and arrays［J］.IEEE Transactions on Signal Processing，2011，59（2）：573－586.

［14］VAIDYANATHAN P P，PAL P.Theory of sparse coprime sensing in multiple dimensions［J］.IEEE Transactions on Signal Processing，2011，59（8）：3592－3608.

［15］WILLETT P K.Modulo conversion method for estimating the direction of arrival［J］.IEEE Transactions on Aerospace and Electronic Systems，2000，36（4）：1391－1396.