本發(fā)明涉及采用陣列天線估計(jì)接收到的信號的到達(dá)方向的技術(shù)領(lǐng)域，尤其涉及采用三維天線陣列的信號到達(dá)方向估計(jì)方法。

背景技術(shù)：

空間信號到達(dá)方向(Direction of Arrival,DOA)估計(jì)是空間譜估計(jì)一個(gè)主要研究方向，被廣泛應(yīng)用在雷達(dá)、聲吶、地震、通信等許多領(lǐng)域。DOA估計(jì)的基本問題就是確定同時(shí)處在空間某一區(qū)域內(nèi)多個(gè)感興趣的信號的空間位置(即各個(gè)信號到達(dá)陣列參考陣元的方向角，簡稱波達(dá)方向)。經(jīng)典的超分辨率DOA估計(jì)算法有多重信號分類算法(MUSIC,Multiple Signal Classification)和基于旋轉(zhuǎn)不變技術(shù)的信號參數(shù)估計(jì)算法(ESPRIT，Estimation of Signal Parameter via Rotational Invitation Techniques)。它們都屬于子空間類算法，其中MUSIC算法是噪聲子空間類算法，ESPRIT算法是信號子空間類算法，以MUSIC算法為代表的算法包括特征矢量法、求根MUSIC法等，以ESPRIT算法為代表的算法包括最小二乘ESPRIT，總體最小二乘ESPRIT等。其中MUSIC算法的中心思想為：利用不同特征值的特征向量之間的正交性將空間劃分為正交的子空間，然后使用這種正交性構(gòu)造陣列空間譜函數(shù)，搜索其極值就可以實(shí)現(xiàn)空間信號電磁波的來向估計(jì)。

傳統(tǒng)的MUSIC算法和ESPRIT算法等高分辨率算法，雖然具有良好的估計(jì)性能，但是由于需要進(jìn)行譜峰搜索或者是對接收信號協(xié)方差矩陣進(jìn)行特征值分解，在應(yīng)用到二維DOA估計(jì)時(shí)具有較大的計(jì)算量，尤其是在陣元數(shù)目較大時(shí)。傳播算子算法在求解信號子空間和噪聲子空間時(shí)僅需要進(jìn)行線性運(yùn)算，因此具有較低的計(jì)算復(fù)雜度。目前，有大量的基于傳播算子的二維DOA估計(jì)算法。但是這些算法具有或者需要進(jìn)行方位角和俯仰角的配對、或者沒有充分利用陣列的結(jié)構(gòu)特點(diǎn)、或者在俯仰角為70°～90°的實(shí)際移動通信俯仰角度范圍內(nèi)存在角度估計(jì)失效問題、或者仍需要進(jìn)行譜峰搜索等局限性。

技術(shù)實(shí)現(xiàn)要素：

為克服現(xiàn)有技術(shù)的不足，本發(fā)明旨在解決傳統(tǒng)的二維DOA估計(jì)中，在俯仰角為70°～90°的實(shí)際移動通信俯仰角度范圍內(nèi)的角度估計(jì)失效問題，以及俯仰角和方位角的配對問題。本發(fā)明采用的技術(shù)方案是，基于三維天線陣列的改進(jìn)傳播算子二維DOA估計(jì)方法，天線陣列為三平行線陣，其中在y軸上有一個(gè)陣元數(shù)目為N+1的均勻線陣Y，在xoy平面內(nèi)和yoz平面內(nèi)分別有一個(gè)和y軸平行的陣元數(shù)目為N的均勻線陣X和Z，子陣的陣元間距，子陣Y和子陣Z以及子陣Y和子陣X之間的距離均為來波信號波長的一半；具體步驟如下：

步驟1：構(gòu)造接收信號矩陣。

將位于坐標(biāo)原點(diǎn)的陣元作為參考陣元，線陣X,Y,Z接收到的信號向量分別為x(t)＝[x₁(t),x₂(t),…,x_N(t)]^T，y(t)＝[y₁(t),y₂(t),…,y_N+1(t)]^T和z(t)＝[z₁(t),z₂(t),…,z_N(t)]^T，其中x_i(t),y_i(t),z_i(t)分別為子陣X,Y,Z的第i個(gè)陣元在t時(shí)刻接收到的信號，構(gòu)造新的接收信號向量w(t)＝[y^T(t),x^T(t),z^T(t)]^T，則對應(yīng)M快拍的接收數(shù)據(jù)矩陣為W＝[w(1),w(2),…,w(M)]；

步驟2：構(gòu)造傳播矩陣

接收信號的自相關(guān)矩陣為將其按如下形式分塊，R＝[R₁,R₂]，其中R₁∈C^(3N+1)×K，R₂∈C^{(3N+1)×(3N+1-K)}，則傳播矩陣定義擴(kuò)展傳播矩陣其中I_K×K為K階單位矩陣，H表示共軛轉(zhuǎn)置；

步驟3：估計(jì)旋轉(zhuǎn)矩陣

將P_c按如下形式分塊，其中P_y∈C^(N+1)×K,P_x∈C^N×K,P_z∈C^N×K，C為復(fù)數(shù)，定義矩陣Ψ_x＝P₁⁺P_x，其中P₁為P_y的前N行，對Ψ_x進(jìn)行特征值分解，則其特征值即為的對角線分量，為Φ_x的估計(jì)值，特征向量矩陣即為A₁的估計(jì)值，其中Φ_x為對應(yīng)子陣X的旋轉(zhuǎn)矩陣，A₁為子陣Y的陣列流型矩陣的前K行；定義兩個(gè)新的矩陣其中P₂為P_y的后N行，P₃為P_x的前N-1行，P₄為P_x的后N-1行，則有其中為Φ_y的估計(jì)值，Φ_y為對應(yīng)子陣Y的旋轉(zhuǎn)矩陣；

同理，定義矩陣Ψ_z＝P₁⁺P_z,對Ψ_z進(jìn)行特征值分解，則其特征值即為的對角線分量，為Φ_z的估計(jì)值，其中Φ_z為對應(yīng)子陣Z的旋轉(zhuǎn)矩陣，特征向量矩陣為即為A₁的估計(jì)值。定義兩個(gè)矩陣其中P₅為P_z的前N-1行，P₆為P_z的后N-1行，則有其中為Φ_y的估計(jì)值；

步驟4：方位角和俯仰角估計(jì)

分別將和按照它們對角線上元素幅角從大到小的順序，對其對角線元素進(jìn)行重新排列，得到新的矩陣和且有其中Π₁和Π₂均為K×K的置換矩陣，令設(shè)分別為的第k個(gè)對角線分量，則方位角和俯仰角的估計(jì)值分別為其中angle表示取幅角運(yùn)算，atan表示取反正切運(yùn)算。

本發(fā)明的特點(diǎn)及有益效果是：

由于本發(fā)明采用基于三維天線陣列的改進(jìn)傳播算子二維估計(jì)算法，因而能夠以較低的計(jì)算復(fù)雜度獲得較好的方位角和俯仰角估計(jì)性能；能夠?qū)崿F(xiàn)方位角和俯仰角估計(jì)的自動配對；在俯仰角為70°～90°的實(shí)際移動通信的俯仰角度范圍內(nèi)不會出現(xiàn)角度模糊問題。

由于本發(fā)明采用的陣列不存在重疊陣元，所以本發(fā)明構(gòu)造的接收信號向量w(t)＝[y^T(t),x^T(t),z^T(t)]^T中不存在冗余數(shù)據(jù)，降低了計(jì)算復(fù)雜度。

本發(fā)明針對三平行線陣，提出了一種解決旋轉(zhuǎn)矩陣之間的配對問題的簡單方法。

附圖說明：

圖1三維天線陣列結(jié)構(gòu)示意圖。

圖2方位角估計(jì)直方圖。

圖3俯仰角估計(jì)直方圖。

圖4不同角度組合估計(jì)聯(lián)合均方誤差圖。

具體實(shí)施方式

針對已有DOA估計(jì)算法存在的問題，本發(fā)明提出了一種基于三維天線陣列的改進(jìn)傳播算子二維DOA估計(jì)算法，該天線陣列為三平行線陣，其中在y軸上有一個(gè)陣元數(shù)目為N+1的均勻線陣Y，在xoy平面內(nèi)和yoz平面內(nèi)分別有一個(gè)和y軸平行的陣元數(shù)目為N的均勻線陣X和Z，子陣的陣元間距，子陣Y和子陣Z以及子陣Y和子陣X之間的距離均為來波信號波長的一半。

本發(fā)明采用的技術(shù)方案：基于三維天線陣列的改進(jìn)傳播算子二維DOA估計(jì)方法，包括以下步驟：

步驟1：構(gòu)造接收信號矩陣。

將位于坐標(biāo)原點(diǎn)的陣元作為參考陣元，線陣X,Y,Z接收到的信號向量分別為x(t)＝[x₁(t),x₂(t),…,x_N(t)]^T，y(t)＝[y₁(t),y₂(t),…,y_N+1(t)]^T和z(t)＝[z₁(t),z₂(t),…,z_N(t)]^T。構(gòu)造新的接收信號向量w(t)＝[y^T(t),x^T(t),z^T(t)]^T，則對應(yīng)M快拍的接收數(shù)據(jù)矩陣為W＝[w(1),w(2),…,w(M)]。

步驟2：構(gòu)造傳播矩陣

接收信號的自相關(guān)矩陣為將其按如下形式分塊，R＝[R₁,R₂]，其中R₁∈C^(3N+1)×K，R₂∈C^{(3N+1)×(3N+1-K)}。則傳播矩陣定義擴(kuò)展傳播矩陣

步驟3：估計(jì)旋轉(zhuǎn)矩陣

將P_c按如下形式分塊，其中P_y∈C^(N+1)×K,P_x∈C^N×K,P_z∈C^N×K。定義矩陣Ψ_x＝P₁⁺P_x，其中P₁為P_y的前N行。對Ψ_x進(jìn)行特征值分解，則其特征值即為的對角線分量，為Φ_x的估計(jì)值，特征向量矩陣即為A₁的估計(jì)值。定義兩個(gè)新的矩陣其中P₂為P_y的后N行，P₃為P_x的前N-1行，P₄為P_x的后N-1行。則有其中為Φ_y的估計(jì)值。

同理定義矩陣Ψ_z＝P₁⁺P_z,對Ψ_z進(jìn)行特征值分解，則其特征值即為的對角線分量，為Φ_z的估計(jì)值，特征向量矩陣為即為A₁的估計(jì)值。定義兩個(gè)矩陣其中P₅為P_z的前N-1行，P₆為P_z的后N-1行。則有其中為Φ_y的估計(jì)值。

步驟4：方位角和俯仰角估計(jì)

分別將和按照它們對角線上元素幅角從大到小的順序，對其對角線元素進(jìn)行重新排列，得到新的矩陣和且有其中Π₁和Π₂均為K×K的置換矩陣。令設(shè)分別為的第k個(gè)對角線分量，則方位角和俯仰角的估計(jì)值分別為其中angle表示取幅角運(yùn)算，atan表示取反正切運(yùn)算。

下面將結(jié)合附圖和實(shí)施例對本發(fā)明做進(jìn)一步的描述：

構(gòu)造如圖1所示的三維天線陣列。假設(shè)空間中有K個(gè)窄帶非相關(guān)信號入射到陣列上，其中第k個(gè)信號的二維波達(dá)方向?yàn)?k＝1,2,…K)，和θ_k分別為來波信號的方位角和俯仰角。

步驟1：構(gòu)造接收信號矩陣。

子陣X,Y,Z在t時(shí)刻接收到的信號向量x(t)＝[x₁(t),x₂(t),…,x_N(t)]^T，y(t)＝[y₁(t),y₂(t),…,y_N+1(t)]^T，z(t)＝[z₁(t),z₂(t),…,z_N(t)]^T可以用式(1)表示。

其中n_x(t)∈C^N×1,n_y(t)∈C^(N+1)×1,n_z(t)∈C^N×1是均值為0，方差為σ²的加性高斯白噪聲，并和s(t)相互獨(dú)立，s(t)為來波信號。為子陣Y的陣列流型矩陣，A_y為A_x的前N行，Φ_x,Φ_z為包含方位角和俯仰角信息的K×K的對角矩陣，且具有如下的表達(dá)形式，

將x(t),y(t),z(t)組合為一個(gè)新的接收信號向量w(t)＝[y^T(t),x^T(t),z^T(t)]^T。則對于M快W＝[w(1),w(2),…,w(M)]。

步驟2：構(gòu)造傳播矩陣

接收信號的自相關(guān)矩陣為將其按如下形式分塊，R＝[R₁,R₂]，其中R₁∈C^(3N+1)×K，R₂∈C^{(3N+1)×(3N+1-K)}。則傳播矩陣定義擴(kuò)展傳播矩陣

步驟3：估計(jì)旋轉(zhuǎn)矩陣

將P_c按如下形式分塊，其中P_y∈C^(N+1)×K,P_x∈C^N×K,P_z∈C^N×K。定義矩陣Ψ_x＝P₁⁺P_x，其中P₁為P_y的前N行。對Ψ_x進(jìn)行特征值分解，則其特征值即為的對角線分量，為Φ_x的估計(jì)值，特征向量矩陣為即為A₁的估計(jì)值。定義兩個(gè)新的矩陣其中P₂為P_y的后N行，P₃為P_x的前N-1行，P₄為P_x的后N-1行。則有為Φ_y的估計(jì)值，其中

同理定義矩陣Ψ_z＝P₁⁺P_z,對Ψ_z進(jìn)行特征值分解，則其特征值即為的對角線分量，為Φ_z的估計(jì)值，特征向量矩陣即為A₁的估計(jì)值。定義兩個(gè)矩陣其中P₅為P_z的前N-1行，P₆為P_z的后N-1行。則有其中為Φ_y的估計(jì)值。

步驟4：方位角和俯仰角估計(jì)

分別將和按照它們對角線上元素幅角從大到小的順序，對其對角線元素進(jìn)行重新排列，得到新的矩陣和且有其中Π₁和Π₂均為K×K的置換矩陣。令設(shè)分別為的第k個(gè)對角線分量，則方位角和俯仰角的估計(jì)值分別為其中angle表示取幅角運(yùn)算，atan表示取反正切運(yùn)算。

結(jié)合上述步驟中的實(shí)施方式，對本發(fā)明的有效性進(jìn)行仿真驗(yàn)證如下：

仿真中取N＝5，即三平行線陣共有16個(gè)陣元，陣列間距d＝0.5λ，其中λ為信號波長，對于每次仿真實(shí)驗(yàn)取快拍數(shù)為200，進(jìn)行M＝500次蒙特卡洛仿真。

仿真實(shí)驗(yàn)1：假設(shè)有K＝2個(gè)等功率非相關(guān)信號入射到天線陣列，其中信噪比SNR＝15dB，信號的方位角和俯仰角為圖2和圖3顯示了方位角估計(jì)直方圖和俯仰角估計(jì)直方圖。從圖中可以看出，本文提出的算法能夠清晰的分辨這兩個(gè)來波信號。

仿真實(shí)驗(yàn)2：假設(shè)有K＝1個(gè)信號入射到天線陣列，信噪比SNR＝15dB，其中信號的方位角和俯仰角均在10°～80°之間以5°的步長變化。圖4為不同角度組合估計(jì)聯(lián)合均方誤差圖。