章鴻運(yùn), 栗蘋, 李國(guó)林, 賈瑞麗

(北京理工大學(xué) 機(jī)電動(dòng)態(tài)控制重點(diǎn)實(shí)驗(yàn)室, 北京 100081)

0 引言

近年來(lái),攜帶軌道角動(dòng)量(OAM)的渦旋電磁波被廣泛應(yīng)用于雷達(dá)[1]、通信[2]和6G無(wú)線網(wǎng)絡(luò)[3]中。電磁波的總角動(dòng)量可以分解為自旋角動(dòng)量和OAM,攜帶不同OAM模式數(shù)的電磁波之間是正交的,這為模分復(fù)用多輸入多輸出(MIMO)系統(tǒng)提供了理論基礎(chǔ)[4]。電磁波的OAM信息在可見光波段具有較為成熟的應(yīng)用。直到2007年渦旋電磁波才首次被引入微波波段中[5]。隨后,使用均勻圓環(huán)陣列產(chǎn)生渦旋電磁波的方法被廣泛采用[6]。

與傳統(tǒng)電磁波不同,渦旋電磁波的波前呈螺旋形分布,提供了方位分辨所需的相位變化。不同模式渦旋電磁波的相位分布不同,實(shí)現(xiàn)了對(duì)角度的分集,還具備一定的抗干擾能力。采用均勻圓陣發(fā)射攜帶不同OAM模式數(shù)的渦旋電磁波,對(duì)接收回波進(jìn)行處理時(shí),目標(biāo)方位信息可以通過模式數(shù)與目標(biāo)方位角之間的對(duì)偶關(guān)系獲得。文獻(xiàn)[7]通過對(duì)回波信號(hào)OAM域做傅里葉變換,獲得目標(biāo)的方位角信息。緊接著空間譜估計(jì)的方法被引入渦旋電磁波探測(cè)中[8-10],使目標(biāo)的角度分辨率得到了提高。由于渦旋電磁波目標(biāo)回波中帶有非線性的貝塞爾函數(shù)項(xiàng),這類方位角估計(jì)方法均對(duì)回波信號(hào)中的貝塞爾函數(shù)進(jìn)行了近似處理,影響了估計(jì)精度。

目前對(duì)均勻圓陣產(chǎn)生渦旋電磁波的分析大多忽略了陣元間的互耦[11],即假定發(fā)射信號(hào)不會(huì)產(chǎn)生二次輻射。此外,根據(jù)奈奎斯特采樣定理,均勻圓陣的陣元間隔小于或等于半波長(zhǎng)是產(chǎn)生渦旋電磁波的必要條件。當(dāng)陣列天線工作在毫米波或者更高頻段時(shí),天線間的距離較短,將導(dǎo)致陣元間的互耦問題。另外,當(dāng)天線為共形陣列時(shí)[12-13],還將出現(xiàn)陰影效應(yīng)和平臺(tái)效應(yīng),陣列性能受到嚴(yán)重影響?；ヱ顚?duì)陣列的影響表現(xiàn)為導(dǎo)向矢量的失配,大多數(shù)高分辨波達(dá)方向估計(jì)(DOA)估計(jì)算法對(duì)陣列流形的誤差比較敏感,將導(dǎo)致角度估計(jì)性能的退化。

互耦導(dǎo)致的誤差可以通過互耦矩陣補(bǔ)償消除。常見的互耦分析方法有開路電壓法[14]、全波法[15]、接收互阻抗法[16]和自校正方法[17-18]等。為解決均勻圓陣中的互耦問題,文獻(xiàn)[14]基于電路模型,將天線陣元上的激勵(lì)電壓分解為球諧模式[19-20]。文獻(xiàn)[16]通過接收互阻抗法得到了與目標(biāo)俯仰角相關(guān)聯(lián)的互阻抗矩陣。文獻(xiàn)[21]結(jié)合陣列S參數(shù)和接收信號(hào),提出了虛擬阻抗互耦校準(zhǔn)方法,該方法與傳統(tǒng)互耦補(bǔ)償矩陣相比具有更小的計(jì)算量。

本文利用均勻圓陣互耦矩陣的循環(huán)對(duì)稱性,通過相位模式構(gòu)造互耦效應(yīng)的等效電路模型,將陣元間的互耦轉(zhuǎn)化至模式域,分析互耦對(duì)渦旋電磁波方位探測(cè)的影響,并采用有源單元方向圖法進(jìn)行仿真。此外,本文利用貝塞爾函數(shù)的遞推關(guān)系,提出一種基于OAM的改進(jìn)傳播算子方法,分析一些常用的典型算法在本文信號(hào)模型下的估計(jì)性能。

1 OAM模型

在無(wú)線電目標(biāo)探測(cè)中,渦旋電磁波可以由均勻圓環(huán)陣列天線產(chǎn)生。如圖1所示,個(gè)陣元均勻分布在半徑為r的圓周上,數(shù)字1、2、n、N表示陣元序號(hào),φ為入射信號(hào)的方位角,φn為第n個(gè)陣元與第1個(gè)陣元的夾角。以圓心為原點(diǎn)建立球坐標(biāo)系,圓陣位于Oxy面上。為了產(chǎn)生模式數(shù)為l的渦旋電磁波,對(duì)每個(gè)陣元施加相位為lφn=2πl(wèi)(n-1)/N(n=1,2,…,N)、幅值為I0的激勵(lì),相鄰兩個(gè)陣元間的相位為Δφ=2πl(wèi)/N。空間中與原點(diǎn)距離為R的任意一點(diǎn)P0(R,θ,φ)的電場(chǎng)E(r)可以表示為

(1)

式中:r為P0點(diǎn)的位置矢量;j為恒定電流密度矢量;ω為信號(hào)角頻率;μ0為真空磁導(dǎo)率;d為偶極子長(zhǎng)度;k為波數(shù),k=2π/λ,λ為波長(zhǎng);rn為第個(gè)陣元的位置矢量,rn=(Rcosφn,Rsinφn,0)。用|r-rn|≈R和|r-rn|≈R-·rn分別作為遠(yuǎn)場(chǎng)幅度和相位的近似[6],其中為r的單位向量。

圖1 均勻圓陣幾何結(jié)構(gòu)Fig.1 UCA array geometry

離散的均勻圓陣可以看作是連續(xù)圓陣在每個(gè)陣元處饋入激勵(lì),即連續(xù)圓陣的采樣。采樣函數(shù)s(φ)可以表示為一個(gè)周期為2π/N的均勻單位脈沖序列:

(2)

式中:δ(·)表示單位階躍函數(shù)。則歸一化電場(chǎng)可表示為

J(l+qN)(krsinθ)ej(l+qN)φ

(3)

式中:Jl為l階第1類貝塞爾函數(shù)。

離散均勻圓陣用于產(chǎn)生渦旋電磁波時(shí),主項(xiàng)由于采樣會(huì)產(chǎn)生周期延拓,式(3)中的求和項(xiàng)包含無(wú)窮多個(gè)高階模式數(shù)。在渦旋電磁波應(yīng)用中只希望產(chǎn)生主模式,而其他周期延拓分量是多余的。根據(jù)貝塞爾函數(shù)的性質(zhì)可知,當(dāng)貝塞爾函數(shù)的階數(shù)大于其輻角krsinθ時(shí)其幅值Jl(krsinθ)≈0,故模式數(shù)l>kr的項(xiàng)可以忽略。另一方面,根據(jù)奈奎斯特采樣定理,N個(gè)陣元的均勻圓陣最多可以產(chǎn)生L=?(N-1)/2」個(gè)渦旋電磁波模式,?·」表示向下取整。因此可以得到kr5時(shí)波束相位不再呈螺旋狀分布,不能產(chǎn)生渦旋電磁波。

2 互耦分析

由第1節(jié)的結(jié)論可知,用于產(chǎn)生渦旋電磁波的均勻圓陣,其陣元間距小于半波長(zhǎng),且電磁波頻率并非是不變常數(shù)。當(dāng)陣列工作在較高頻段時(shí),產(chǎn)生渦旋電磁波不僅受空間采樣影響,陣元間還存在互耦。

在陣列天線分析中,互耦效應(yīng)常用阻抗矩陣Z描述。考慮所有陣元之間都存在互耦,即互耦自由度為N,則

(4)

式中:Zmn為第m個(gè)和第n個(gè)陣元間的互阻抗,

Zmn=Vm/In

(5)

Vm(m=1,…,N)為第m個(gè)陣元的開路電壓,In(n=1,…,N)為當(dāng)其他陣元開路、第n個(gè)陣元的短路電流。如圖3所示,發(fā)射信號(hào)經(jīng)相位加權(quán)后饋給各個(gè)天線,發(fā)射阻抗ZG1,…,ZGN可以表示為一個(gè)對(duì)角矩陣ZG,則互耦矩陣C為

C=(ZG+Z)-1Z

(6)

圖3中,Port1、Port2、…、PortN分別為天線的等效端口。將上述互耦模型代入渦旋電磁波信號(hào)中,可得

E(θ,φ)=ICa(θ,φ)

(7)

式中:I為第l個(gè)OAM模式的激勵(lì)矢量,

I=I0[ejlφ1,ejlφ2,…,ejlφn]

(8)

a為理想情況下的導(dǎo)向矢量,

(9)

通過以上分析可知,陣列天線的性質(zhì)可以完全由阻抗矩陣Z和a(θ,φ)表征。

圖2 理想情況12陣元不同OAM模式相位分布圖Fig.2 Phase distribution of different OAM modes of 12 elements in ideal conditions

對(duì)于每一個(gè)OAM模式,可以通過對(duì)互耦矩陣做離散傅里葉變換,將各個(gè)陣元間的互耦分解成相位序列c,其中第l個(gè)元素為

(10)

式中:Cmn為第m個(gè)和第n個(gè)陣元間的互耦。當(dāng)圓陣中的陣元相同且各向同性時(shí),任意陣元的阻抗都相同,互耦矩陣C具有循環(huán)對(duì)稱性,式(10)可以簡(jiǎn)化為對(duì)互耦矩陣的第1行做離散傅里葉變換(DFT):

(11)

圖4將均勻圓陣中的互耦表示成OAM模式,其中等效電流,等效電阻和等效電壓分別展開成自變量為l的序列。圖4中,Voc為開路電壓源,zl為等效負(fù)載阻抗。由此可得互耦情況下空間任意一點(diǎn)的歸一化電場(chǎng)為

E(θ,φ)=cljlJl(krsinθ)ejlφ

(12)

3 方位估計(jì)

3.1 渦旋電磁波接收信號(hào)模型

基于多發(fā)單收的渦旋電磁波目標(biāo)探測(cè)原理如圖5所示。振蕩器產(chǎn)生的載頻信號(hào)經(jīng)波形調(diào)制器調(diào)制后,經(jīng)由功率放大器進(jìn)行放大,放大后的射頻信號(hào)經(jīng)過功分器分為功率相同的N路信號(hào),接著移相器給每路信號(hào)附加Δφ=2πl(wèi)/N的相移,再通過均勻圓陣依次發(fā)射[-(2L+1),2L+1]OAM模式數(shù)的渦旋電磁波。接收天線位于陣列圓心,回波信號(hào)進(jìn)入混頻器與載頻信號(hào)進(jìn)行混頻得到中頻信號(hào),經(jīng)放大濾波后在OAM域進(jìn)行采樣,最后通過信號(hào)處理得到目標(biāo)方位角信息。

假設(shè)空間中存在P個(gè)目標(biāo),回波信號(hào)和噪聲均為平穩(wěn)、零均值、不相關(guān)的隨機(jī)過程。則接收回波信號(hào)可以表示為

(13)

圖3 陣列互耦饋電網(wǎng)絡(luò)圖Fig.3 Schematic of the array mutual coupling feeding network

圖4 均勻圓陣互耦OAM模式等效電路圖Fig.4 OAM mode circuit model of the UCA including mutual coupling

圖5 渦旋電磁波探測(cè)原理圖Fig.5 Schematic diagram of vortex electromagnetic wave detection

式中:sp(t)為第p個(gè)空間信號(hào)。寫成矢量形式為

X(t)=ClBS(t)+N(t)

(14)

式中:Cl=diag(c);B為模式空間流形矩陣,B=[b1,b2,…,bp,…,bP],

bp=jlJl(krsinθp)ejlφp,l={-L,…,0,…,L}T

(15)

S(t)為信號(hào)矢量;N(t)為功率為σ2的高斯白噪聲。

3.2 MUSIC算法

接收數(shù)據(jù)的協(xié)方差矩陣為

(16)

式中:RS為信號(hào)協(xié)方差矩陣,RS=E[S(t)SH(t)]。對(duì)R進(jìn)行特征值分解,令US為信號(hào)子空間,則譜估計(jì)[22]為

(17)

式中:W=diag(jlJl(krsinθ));v(φ)=ejlφ。方位角可以通過對(duì)P個(gè)譜峰的搜索得到。圖6為MUSIC算法的方位角估計(jì)譜,其中陣元數(shù)N=12,半徑r=λ,目標(biāo)方位角分別為50°和80°。由圖6可以看出,陣列互耦對(duì)目標(biāo)方位角估計(jì)產(chǎn)生了嚴(yán)重影響,通過互耦校準(zhǔn)可明顯改善陣列估計(jì)性能。

3.3 OAM-PM算法

MUSIC算法需要對(duì)協(xié)方差矩陣進(jìn)行特征值分解和譜峰搜索,計(jì)算量較大。傳播算子算法利用線性運(yùn)算替代了特征值分解,降低了計(jì)算復(fù)雜度?，F(xiàn)有的傳播算子算法在對(duì)導(dǎo)向矢量矩陣分塊時(shí),通常將非奇異矩陣用導(dǎo)向矢量的前P行來(lái)定義[23]。然而在渦旋電磁波探測(cè)中,由式(13)可知不同模式數(shù)的回波信號(hào)幅值受貝塞爾函數(shù)調(diào)制,模式數(shù)較低時(shí)對(duì)應(yīng)的幅值較大。因此本文為了提高方位估計(jì)的穩(wěn)健性,對(duì)傳播算子矩陣進(jìn)行了重新定義,提出一種基于OAM的傳播算子方位估計(jì)方法。

將矩陣B分塊:

(18)

(19)

式中:P為傳播算子矩陣,P∈CP×(2L+1-P)。定義矩陣:

(20)

式中:P1為P的前?(2L+1-P)/2」列;P2為P的后「(2L+1-P)/2?列。由式(19)和式(20)可得

QBy=B

(21)

矩陣Q∈C(2L+1)×P張成的空間與矩陣B張成的模式空間是同一個(gè)空間。

在噪聲和互耦情況下,P可以從接收信號(hào)中得到。計(jì)算協(xié)方差矩陣

(22)

式中:t為時(shí)間變量;K為快拍數(shù)。將分塊

(23)

式中:R1∈C(2L+1)×?(2L+1-P)/2」;R2∈C(2L+1)×P;R3∈C(2L+1)×「(2L+1-P)/2?。根據(jù)式(18)和式(19)可得

(24)

通過最小二乘法可以得到P的估計(jì):

(25)

由于傳播矩陣中含有貝塞爾函數(shù),可以利用式(26)的遞推關(guān)系[24]:

(2l/krsinθ)Jl(krsinθ)=Jl-1(krsinθ)+Jl+1(krsinθ)

(26)

取Q和B的前2L-1行、第2～2L行和后2L-1行,分別定義為Q1、Q2、Q3和B1、B2、B3。將分塊矩陣分別代入式(21),可得

Q1Bx=B1,Q2Bx=B2,Q3Bx=B3

(27)

根據(jù)貝塞爾函數(shù)的性質(zhì)可得B1、B2、B3具有如下關(guān)系:

DB2=B1Φ+B3Φ*

(28)

式中:D=(2/krsinθ)diag{-(L-1),…,0,…,(L-1)};Φ=diag(ej(φ1+π/2),ej(φ2+π/2),…,ej(φP+π/2));Φ*為Φ的共軛轉(zhuǎn)置。同理可得

DQ2=Q1Ψ+Q3Ψ*

(29)

(30)

可得

(31)

式中:?表示Moore-Penrose偽逆;Ψ*為Ψ的共軛轉(zhuǎn)置。對(duì)Ψ特征值分解即可獲得目標(biāo)的方位角:

φp=arg (λp)-π/2

(32)

經(jīng)過以上分析,本文所提解耦合OAM方位估計(jì)算法步驟如下:

2) 將協(xié)方差矩陣分塊,由式(25)得到傳播算子矩陣P,并通過式(20)得到矩陣Q。

3) 將矩陣Q劃分為Q1、Q2、Q3,分別為Q的前2L-1行、第2～2L行和后2L-1行。

4 仿真分析

4.1 陣列互耦仿真

互耦矩陣通過有源單元方向圖法得到[25]。天線的有源單元方向圖是指陣列中只激勵(lì)一個(gè)單元,同時(shí)其他所有單元終端都接匹配負(fù)載時(shí)的輻射方向圖。當(dāng)陣列天線單元間存在互耦時(shí),陣列環(huán)境中的單元輻射方向圖與孤立單元方向圖不同。有源單元方向圖法是將孤立單元方向圖和有源單元方向圖通過互耦矩陣建立聯(lián)系:

(33)

通過最小二乘法求解互耦矩陣:

(34)

陣列參數(shù)由表1所示。

表1 FEKO仿真參數(shù)Table 1 FEKO simulation parameters

通過FEKO仿真得到耦合相位序列為

c=[0.743 1+0.229 0i,0.980 4- 0.041 3i,0.894 1-0.082 5i 1.242 8+0.120 0i,1.280 1-0.261 7i, 1.245 7+0.121 3i 1.280 1-0.261 7i,1.242 8+0.120 0i, 0.894 1-0.082 5i 0.980 4-0.041 3i,0.743 1+0.229 0i]

(35)

4.2 計(jì)算復(fù)雜度

本文算法涉及協(xié)方差矩陣計(jì)算、矩陣特征值分解、矩陣偽逆和最小二乘法等運(yùn)算。假設(shè)快拍數(shù)為K,譜峰搜索范圍為[0°,360°),譜峰搜索步長(zhǎng)為1°,各算法計(jì)算復(fù)雜度如表2所示。

表2 不同算法的計(jì)算復(fù)雜度Table 2 Computational complexity of different algorithms

圖7為快拍數(shù)K=100、最大模式數(shù)L變化范圍為5～100時(shí)算法計(jì)算出的復(fù)雜度曲線。

圖7 不同模式數(shù)的計(jì)算復(fù)雜度Fig.7 Computational complexity of different mode numbers

圖8為模式數(shù)L=50、快拍數(shù)K變化范圍為50～1 000時(shí)算法計(jì)算出的復(fù)雜度曲線。

從圖7和圖8中可以看出,本文算法將MUSIC算法和ESPRIT算法中2L+1維的矩陣特征值分解改進(jìn)為求P維矩陣的Moore-Penrose逆,同時(shí)還避免了譜峰搜索,降低了計(jì)算復(fù)雜度。

4.3 性能分析

為驗(yàn)證所提算法在不同陣列配置下的方位探測(cè)性能,在陣元數(shù)N=12、半徑r=λ的均勻圓陣上進(jìn)行了200次獨(dú)立的蒙特卡洛仿真。其他仿真參數(shù)設(shè)置如下:載頻f=24 GHz,總模式數(shù)2L+1=11,目標(biāo)個(gè)數(shù)P=2,方位角分別為30°和80°。

無(wú)條件信號(hào)模型的克拉美羅界[26]可以表示為

(36)

定義均方根誤差為

(37)

圖9為快拍數(shù)為1 000、信噪比變化范圍為0～30 dB時(shí),不同算法計(jì)算出的均方根誤差隨信噪比變化的曲線。

圖10為信噪比為10 dB、快拍數(shù)變化范圍為50～1 000時(shí),不同算法均方根誤差隨快拍數(shù)變化的曲線。

圖10 均方根誤差與快拍數(shù)關(guān)系圖Fig.10 Dependence of root mean square error on snapshots

由圖10可以看出,本文所提算法相較于傳統(tǒng)傳播算子算法估計(jì)誤差明顯降低,性能與MUSIC算法和ESPRIT算法近似,且具備更低的計(jì)算復(fù)雜度。

5 結(jié)論

本文提出一種解耦合的OAM目標(biāo)方位探測(cè)新算法。首先建立OAM互耦模型,使用FEKO軟件對(duì)12陣元的均勻圓環(huán)陣列偶極子天線仿真得到互耦矩陣,然后通過DFT將其轉(zhuǎn)化到模式域,最后通過OAM-PM算法進(jìn)行目標(biāo)方位估計(jì)。實(shí)驗(yàn)結(jié)果表明,所提算法具有良好的解耦合性能,對(duì)目標(biāo)方位角具備精確探測(cè)能力的同時(shí)擁有更低的計(jì)算復(fù)雜度。