張狀和 韓東 劉德亮

摘? 要： 針對柱面共形陣的二維波達(dá)方向（2D DOA）估計(jì)問題，提出一種基于迭代自適應(yīng)方法的柱面共形陣2D DOA估計(jì)算法。首先，通過歐拉旋轉(zhuǎn)變換公式將柱面共形陣全局直角坐標(biāo)系變換為單元局部直角坐標(biāo)系，建立柱面共形天線陣列流型;然后，將感興趣區(qū)域內(nèi)信號源所有可能存在的位置通過劃分二維網(wǎng)格的形式來表示，利用迭代自適應(yīng)算法（IAA）估計(jì)出每個(gè)潛在位置對應(yīng)的信號源能量;最后，繪制能量譜圖，譜峰對應(yīng)的方位角和俯仰角即為二維波達(dá)方向的估計(jì)值。理論分析與仿真實(shí)驗(yàn)證明，所提算法不需要額外的參數(shù)匹配，簡化了估算步驟，在短快拍的情況下可以實(shí)現(xiàn)柱面共形陣的2D DOA估計(jì)。

關(guān)鍵詞： 柱面共形陣; 2D DOA估計(jì); 坐標(biāo)系變換; 信號源能量估計(jì); 迭代自適應(yīng)算法; 能量譜繪制; 理論分析

中圖分類號： TN911.7?34? ? ? ? ? ? ? ? ? ? ? ?文獻(xiàn)標(biāo)識碼： A? ? ? ? ? ? ? ? ? ? ? ? ? 文章編號： 1004?373X（2020）11?0006?04

Cylindrical conformal array′s 2D DOA estimation based on iterative adaptive approach

ZHANG Zhuanghe， HAN Dong， LIU Deliang

（Department of Missile Engineering， Army Engineering University of PLA （Shijiazhuang Campus）， Shijiazhuang 050003， China）

Abstract： An iterative adaptive approach （IAA） based 2D DOA （two dimensional direction of arrival） estimation algorithm for cylindrical conformal array is proposed to deal with the problem of 2D DOA estimation of cylindrical conformal array. The global rectangular coordinate system of the cylindrical conformal array is transformed into the unit local rectangular coordinate system by Euler rotation transformation formula for the establishment of array manifold of the cylindrical conformal antenna array. All possible locations of signal source in the ROI （region of interest） are represented by dividing 2D grids， and the IAA is used to estimate the signal source energy corresponding to each potential location. An energy spectrum map is drawn， and the azimuth angle and pitch angle corresponding to the spectrum peak are the estimated values of the 2D DOA. Theoretical analysis and simulation experiments prove that it is unnecessary to match additional parameters in the application of the proposed algorithm， so the steps of estimation is simplified and 2D DOA estimation of cylindrical conformal array can be realized under the condition of short snapshot.

Keywords： cylindrical conformal array; 2D DOA estimation; coordinate system transformation; signal source energy estimation; iterative adaptive approach （IAA）; energy spectrum drawing; theoretical analysis

0? 引? 言

共形陣天線是一種和載體外形保持一致的天線[1]，其能夠和物體結(jié)構(gòu)相融合，可以滿足飛機(jī)、導(dǎo)彈、衛(wèi)星等飛行器的隱身和空氣動力學(xué)需求，有效減少了荷載的大小和重量。相對于傳統(tǒng)的線陣和平面陣，共形陣天線還具有覆蓋范圍大、易偽裝、抗干擾能力較強(qiáng)等特點(diǎn)。近些年來，學(xué)者們圍繞著共形陣列展開了多方面的研究[2?7]，其中，由于波達(dá)方向（Direction of Arrival，DOA）估計(jì)技術(shù)在通信、雷達(dá)、聲吶等領(lǐng)域的重要作用，因此如何實(shí)現(xiàn)共形陣列DOA估計(jì)受到廣泛關(guān)注。

文獻(xiàn)[2]研究了多重信號分類（Multiple Signal Classification，MUSIC）算法在共形天線中的應(yīng)用。文獻(xiàn)[3]通過虛擬陣列變換技術(shù)，應(yīng)用旋轉(zhuǎn)不變性技術(shù)估計(jì)信號參數(shù)（Estimation of Signal Parameters Via Rotational Invariance Techniques，ESPRIT）和MUSIC算法，克服了共形載體的遮蔽效應(yīng)，實(shí)現(xiàn)了基于共形陣的2D DOA（Two Dimensional Direction of Arrival，2D DOA）估計(jì)。但是以上算法忽略了載體曲率對方向圖綜合的影響，簡化了陣列流型的數(shù)據(jù)模型。文獻(xiàn)[4?5]通過運(yùn)用歐拉旋轉(zhuǎn)變換公式，提出了共形天線陣列流型的通用建模方法。文獻(xiàn)[6]通過合理設(shè)置陣元，并基于ESPRIT，解決了柱面共形陣的盲極化2D DOA估計(jì)。文獻(xiàn)[7]則引入了虛擬期望信號，利用最大信干噪比準(zhǔn)則下的最優(yōu)權(quán)矢量對接收信號進(jìn)行加權(quán)處理，實(shí)現(xiàn)了共形陣DOA估計(jì)。上述文獻(xiàn)均將經(jīng)典的子空間類算法應(yīng)用到了共形陣列DOA估計(jì)問題上，但是這些算法都需要大量的快拍數(shù)據(jù)支撐，以保證對樣本協(xié)方差矩陣進(jìn)行特征分解時(shí)，信號子空間與噪聲子空間不發(fā)生混疊，同時(shí)獲得高精度、高分辨率的估計(jì)結(jié)果。但是在實(shí)際的工程應(yīng)用中，信號經(jīng)常處于不能長時(shí)間穩(wěn)定或快速時(shí)變的環(huán)境中，如高速目標(biāo)追蹤、水下信號處理、城市無線通信等，采集大量的快拍數(shù)據(jù)會導(dǎo)致處理時(shí)間過長，增大與真實(shí)樣本的誤差或降低運(yùn)算速度等問題。而近些年來，隨著稀疏信號處理理論的興起，尤其以Petre Stoica等提出的振幅和相位估計(jì)的迭代自適應(yīng)方法（Iterative Adaptive Approach for Amplitude and Phase Estimation， IAA?APES）[8?9]、基于協(xié)方差稀疏迭代的譜估計(jì)方法（Sparse Iterative Covariance?based Estimation Method， SPICE）[10]為代表，使得短快拍DOA估計(jì)的實(shí)現(xiàn)成為可能。

為了能在短快拍條件下實(shí)現(xiàn)共形陣列天線的高精度2D DOA估計(jì)，本文提出了一種基于迭代自適應(yīng)方法的柱面共形陣2D DOA估計(jì)算法。首先，利用歐拉旋轉(zhuǎn)變換公式，實(shí)現(xiàn)柱面共形陣天線全局直角坐標(biāo)系到單元局部直角坐標(biāo)系的變換，建立柱面共形天線陣列流型;然后，基于加權(quán)最小二乘法估計(jì)出入射信號源的能量，通過循環(huán)迭代對估計(jì)結(jié)果進(jìn)行更新直至收斂;最后，繪制出三維能量譜圖，譜峰對應(yīng)方位角和俯仰角即為入射信源的2D DOA。仿真實(shí)驗(yàn)表明，在只有2個(gè)快拍的情況下，所提算法經(jīng)過適量迭代就能實(shí)現(xiàn)對入射信號源2D DOA的高精度、高分辨率估計(jì)。

下面對文中的符號和算子進(jìn)行說明：[[?]*]表示復(fù)數(shù)共軛;[[?]T]表示轉(zhuǎn)置;[[?]H]表示共軛轉(zhuǎn)置;[?]表示克羅內(nèi)克積;[x]表示[x]的估計(jì)值;[diag（? ）]表示對角化;[E[?]]表示期望值;[?]表示Euclidean范數(shù)。

1? 柱面共形陣數(shù)據(jù)模型

柱面共形陣列天線結(jié)構(gòu)如圖1所示。[M×N] 個(gè)陣元均勻地分布在載體表面，定義[n]為由下至上的陣元序號，[m]為每一圓環(huán)上逆時(shí)針方向的陣元序號，[R]表示圓柱半徑，[H]表示圓柱的高，[d]為上下兩陣元間的距離，[β] 為相鄰兩陣元間的夾角，[?]為圓環(huán)末陣元與[x]軸的夾角;并定義坐標(biāo)系[[x，y，z]]表示陣列全局直角坐標(biāo)系，[[x，y，z]] 表示陣元局部直角坐標(biāo)系，[[r，θ，φ]]表示陣列全局極坐標(biāo)，[[r，θ，φ]]表示陣元局部極坐標(biāo)。

給出陣元坐標(biāo)以及歐拉旋轉(zhuǎn)變換角公式[4]：

[xnm=Rcos[（m-1）β-?]] （1）

[ynm=Rsin[（m-1）β-?]] （2）

[znm=nd-（N+1）d/2] （3）

[Dnm=π-?+（M-1）βEnm=-π2Fnm=0] （4）

假設(shè)有[K]個(gè)遠(yuǎn)場窄帶信號源入射到圖1所示柱面共形陣上，[θi]，[?i]分別表示第[i]個(gè)入射信號源的方位角與俯仰角，結(jié)合文獻(xiàn)[4?5]，柱面共形陣的導(dǎo)向矢量模型可以表示為：

[a（θ，?）=[r1e-j2πp1?uλ，r2e-j2πp2?uλ，…，rNMe-j2πpnm?uλ]T] （5）

[ri=（g2iθ+g2i?）12（k2iθ+k2i?）12cos θigk=giplcos θigk=gi?pl=giθkθ+gi?k?] （6）

式中：[pnm=[xnm，ynm，znm]]（[n=1，2，…，N]，[m=1，2，…，M]）表示第[nm] 個(gè)陣元的坐標(biāo)矢量;如圖2所示，[u=[cos θsin ?，sin θsin ?，cos ?]T]定義為傳播矢量;[ri]表示第[i]個(gè)陣元對單位信號的響應(yīng);[gi]為陣元方向圖表達(dá)式;[giθ，gi?]分別由信源構(gòu)成的矢量基[uθ，u?]所構(gòu)成;[θigk]為[gi]和[pl]向量之間的夾角。

考慮到入射信號源均相互獨(dú)立，噪聲為加性高斯白噪聲，因此快拍數(shù)據(jù)模型可以表示為：

[X（n）=AS（n）+N（n）] （7）

[S（n）=[s1，s2，…，sn]T] （8）

[N（n）=[n1，n2，…，nn]T] （9）

[A（θ，?）=[a（θ1，?1），a（θ2，?2），…，a（θNM，?NM）]] （10）

式中：[A（θ，?）]表示導(dǎo)向矩陣，由導(dǎo)向矢量構(gòu)成;[S]表示信號矢量;[N]表示加性高斯白噪聲矢量。通常情況下，實(shí)際感興趣的目標(biāo)信號源數(shù)大于掃描區(qū)域內(nèi)信號源的個(gè)數(shù)，導(dǎo)致[S（n）]中只有一小部分的數(shù)值是非零的，因此稀疏類算法適用于DOA估計(jì)。

2? 柱面共形陣2D DOA估計(jì)算法

首先，構(gòu)造[K×K]的能量矩陣[P]，其對角線元素包含了掃描網(wǎng)格上每個(gè)角度的信號能量。能量矩陣的結(jié)構(gòu)如下：

[P=diag{p1，p2，…，pi}] （11）

[pi=1Nn=1Nsi（n）2，? ? i=1，2，…，K]? （12）

定義噪聲和干擾的協(xié)方差矩陣如下：

[Q（θi，?i）=R-pia（θi，?i）aH（θi，?i）]? ?（13）

式中[R=A（θ，?）PAH（θ，?）]。

加權(quán)最小二乘法的代價(jià)函數(shù)由下式給出：

[n=1Nx（n）-si（n）a（θi，?i）2W] （14）

式中：[y2W=yHWy];[si（n）]是第[i]個(gè)遠(yuǎn)場信號源在快拍數(shù)[n]處的復(fù)幅值;[Q-1（θi，?i）]是噪聲的協(xié)方差矩陣，由文獻(xiàn)[7]可知，當(dāng)[W=Q-1（θi，?i）]時(shí)[si（n）]的值取最小，因此，對式（14）求關(guān)于[si（n）]的最小化，表達(dá)形式為：

[si（n）=aH（θi，?i）Q-1（θi，?i）x（n）aH（θi，?i）Q-1（θi，?i）a（θi，?i）， n=1，2，…，N] （15）

式中：[Q-1（θi，?i）]可由[R-1]代替，式（15）可以精簡為：

[si（n）=aH（θi，?i）R-1x（n）aH（θi，?i）R-1a（θi，?i），? ?n=1，2，…，N] （16）

式（16）得到了[（θi，?i）]處對應(yīng)信源復(fù)幅值的估計(jì)值。將[si（n）]代入式（12）即可計(jì)算第[i]個(gè)掃描點(diǎn)的能量。

算法流程總結(jié)如下：

初始化：

[s（0）i（n）=aH（θi，?i）x（n）a（θi，?i）2，? ? i=1，2，…，K] （17）

[p（0）i=1Nn=1Ns（0）i（n）2，? ? i=1，2，…，K] （18）

[P（0）=diag{p（0）1，p（0）2，…，p（0）K}] （19）

迭代過程如下：

1） 構(gòu)造能量矩陣[p（l）i=1Nn=1Ns（l-1）i（n）2]，[i=1，][2，…，K];[p（l）i]是第[i]個(gè)網(wǎng)格點(diǎn)對應(yīng)的信號源在第[l]次迭代的能量估計(jì)值。

2） 計(jì)算信號協(xié)方差矩陣[R=A（θ，?）P（l）AH（θ，?）]。

3） 更新[s（l）i（n）=aH（θi，?i）R-1x（n）aH（θi，?i）R-1a（θi，?i）]。

4） 當(dāng)?shù)諗?，則停止迭代，否則返回式（1），根據(jù)經(jīng)驗(yàn)，一般迭代10～15次就會得到較好的估計(jì)值。

計(jì)算能量值：迭代結(jié)束后，以最后一次迭代得到的[p（end）i]作為第[i]個(gè)掃描點(diǎn)對應(yīng)的譜峰函數(shù)繪制出三維譜峰圖，譜峰對應(yīng)的方位角、俯仰角即為真實(shí)信源的二維波達(dá)角度。

3? 仿真實(shí)驗(yàn)與分析

為驗(yàn)證本文所提方法的有效性，對比2D?MUSIC算法展開實(shí)驗(yàn)。為了不失一般性，設(shè)置如圖1所示的接收陣列，[M=5]，[N=10]，相鄰兩圓環(huán)的間距為入射信號波長的一半，即[d=λ/2];噪聲為加性高斯白噪聲，與入射信號源相互獨(dú)立;信號源個(gè)數(shù)[K=2]，且相互獨(dú)立;[θ1=55°]，[?1=30°]，[θ2=60°]，[?2=35°];假設(shè)在DOA估計(jì)之前已對感興趣的區(qū)域有一定認(rèn)知，所以提前設(shè)定方位角度掃描范圍為[[50°，65°]]，俯仰角度掃描范圍為[[25°，40°]]，掃描間隔為[0.5°];迭代次數(shù)預(yù)設(shè)[l=10]。單元方向圖設(shè)置為[4]：

當(dāng)[0≤θ≤π2]時(shí)：

[gθ（θ，?）=J2πdλsin θ-J0πdλsin θ×（cos ?-jsin ?）]? （20）

[g?（θ，?）=J2πdλsin θ+J0πdλsin θ×cos θ（sin ?-jcos ?）] （21）

其他時(shí)，

[gθ（θ，?）=g?（θ，?）=0]

式中[J0，J2]分別是零階和二階第一類貝塞爾函數(shù)。

3.1? 仿真實(shí)驗(yàn)1

設(shè)置信噪比（Signal to Noise Ratio，SNR）為15 dB，快拍數(shù)[N1=2，N2=200]，能量譜分別如圖3，圖4所示。

圖3，圖4仿真結(jié)果表明，當(dāng)快拍數(shù)足夠大時(shí)，兩種算法都能成功地把三個(gè)近場信源分辨出來;而當(dāng)快拍數(shù)減小到[N1=2] 時(shí)，2D?MUSIC算法已經(jīng)失效，真實(shí)信號源的譜峰混疊在偽峰之中，已無法分辨出真實(shí)入射角度，但所提算法仍然可以分辨出兩個(gè)銳利的譜峰，且旁瓣水平較低，方位角和俯仰角的估計(jì)偏差也很小。

3.2? 仿真實(shí)驗(yàn)2

定義實(shí)驗(yàn)成功須滿足下式：

[（θ）-θ≤3°，? ? （?）-?≤3°] （22）

圖5設(shè)置快拍數(shù)[N2=200]，信噪比從-15 dB變化到15 dB。圖6設(shè)置信噪比為10 dB，快拍數(shù)從10變化到200，其他實(shí)驗(yàn)條件不變。200次蒙特卡洛實(shí)驗(yàn)觀察分辨概率和均方根誤差隨快拍數(shù)的變化曲線。

從圖5的實(shí)驗(yàn)結(jié)果可以看出，所提算法的估計(jì)均方根誤差要小于2D?MUSIC算法，說明所提算法的估計(jì)精度更高。從圖6的實(shí)驗(yàn)結(jié)果可以看出，當(dāng)快拍數(shù)大于70時(shí)，所提算法的分辨概率趨于100%，而2D?MUSIC算法的分辨概率依舊很低。

4? 結(jié)? 語

針對柱面共形陣的2D DOA估計(jì)問題，提出一種基于IAA的估計(jì)算法。本文利用歐拉旋轉(zhuǎn)變換公式，將柱面共形陣天線全局直角坐標(biāo)系變換到單元局部直角坐標(biāo)系，討論了稀疏類算法在DOA估計(jì)應(yīng)用上的原理，使用IAA算法解決了短快拍條件下的柱面共形陣的2D DOA估計(jì)問題，計(jì)算機(jī)仿真驗(yàn)證了所提算法的有效性。

參考文獻(xiàn)

[1] JOSEFSSO L， PERSSON P. Conformal array antenna theory and design [M]. New Jersey： John Wiley & Sons， Inc.， 2006.

[2] 孫翠珍，毛昕蓉.多重信號分類算法在共形陣列天線中的應(yīng)用[J].計(jì)算機(jī)仿真，2014，31（11）：171?174.

[3] 王瑞革，王法棟.基于虛擬陣列變換的共形陣列信號DOA估計(jì)[J].雷達(dá)科學(xué)與技術(shù)，2016，14（4）：448?452.

[4] 王布宏，郭英，王永良，等.共形天線陣列流形的建模方法[J].電子學(xué)報(bào)，2009，37（3）：481?484.

[5] 李龍軍，王布宏，夏春和.稀疏共形陣列天線方向圖綜合[J].電子學(xué)報(bào)，2017，45（1）：104?111

[6] 張羚，郭英，齊子森，等.柱面共形陣列天線盲極化2D DOA估計(jì)[J].空軍工程大學(xué)學(xué)報(bào)（自然科學(xué)版），2016，17（3）：78?84.

[7] 楊群，曹祥玉，高軍，等.一種未知信源數(shù)的共形陣DOA算法[J].西安電子科技大學(xué)學(xué)報(bào)（自然科學(xué)版），2015，42（3）：122?128.

[8] YARDIBI T， LI J， STOICA P， et al. Source localization and sensing： a nonparametric iterative adaptive approach based on weighted least squares [J]. IEEE transactions on aerospace and electronic system， 2010， 46（1）： 425?443.

[9] STOICA P， LI J， HE H. Spectral analysis of nonuniformly sampled data： a new approach versus the periodogram [J]. IEEE transactions on signal processing， 2009， 57（3）： 843?858.

[10] STOICA P， BABU P， LI J. New method of sparse parameter estimation in separable models and its use for spectral analysis of irregularly sampled data [J]. IEEE transactions on signal processing， 2011， 59（1）： 35?47.