本發(fā)明屬于空間譜估計(jì)領(lǐng)域，具體涉及一種高采樣1比特量化情況下的信號(hào)到達(dá)角高精度估計(jì)方法。

背景技術(shù)：

空間譜估計(jì)是利用在空間中，以不同排列方式組合在一起的接收陣列，對(duì)信源信息進(jìn)行接收并處理，從而獲得信源參數(shù)的技術(shù)?？臻g譜估計(jì)對(duì)信源空域參數(shù)，如信源個(gè)數(shù)，信號(hào)波達(dá)方向(DOA)及俯仰角的優(yōu)越估計(jì)性能，已經(jīng)使其廣泛應(yīng)用在通信，雷達(dá)，聲吶及生物醫(yī)學(xué)等領(lǐng)域。常見的空間陣列排列方式有，均勻線陣，均勻圓陣，平面陣。其中，均勻線陣在實(shí)際分析中應(yīng)用廣泛。本發(fā)明后續(xù)采用均勻線陣進(jìn)行實(shí)例的仿真實(shí)驗(yàn)。

對(duì)于信號(hào)波達(dá)方向(DOA)的估計(jì)問題，有很多學(xué)者從不同方面提出了空間譜估計(jì)的方法，有ARMA譜分析法，最大似然法，熵譜分析法和特征分解法等，其中基于特征分解的子空間算法主要有MUSIC，ESPRIT，WSF等，但是這些算法在有些信號(hào)模型下會(huì)失效，如信源是相干的，并且在實(shí)際環(huán)境中信噪比低，采樣數(shù)不夠多，獲取的自相關(guān)矩陣不準(zhǔn)確時(shí)，這些算法的估計(jì)性能也會(huì)大大下降。近年來，隨著稀疏理論的不斷完善，利用稀疏信號(hào)表達(dá)的方法對(duì)未知信號(hào)波達(dá)方向的估計(jì)已經(jīng)取得了進(jìn)步，避免了自相關(guān)矩陣的運(yùn)算，也不受信源相干性的影響，同時(shí)也提高了估計(jì)性能。

在實(shí)際應(yīng)用中，陣列接收機(jī)接收到的數(shù)據(jù)需要被量化。1比特采樣量化就是利用1位二進(jìn)制數(shù)據(jù)來表示觀測(cè)信號(hào)中的符號(hào)信息，在工程實(shí)現(xiàn)中，也不需要復(fù)雜的電路結(jié)構(gòu)，只需要一個(gè)比較器就可以完成，大大降低了工程實(shí)現(xiàn)成本。

因此，采用1比特采樣量化并稀疏重構(gòu)原始信號(hào)信息的方法比傳統(tǒng)意義上的DOA估計(jì)算法，估計(jì)精確要好，采樣要求低，工程實(shí)現(xiàn)簡(jiǎn)單，快速。

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

本發(fā)明的目的是提供一種高采樣1比特量化情況下的信號(hào)到達(dá)角高精度估計(jì)方法，解決了在低信噪比下，采樣量少，自相關(guān)矩陣獲取不準(zhǔn)確的問題。

本發(fā)明所采用的技術(shù)方案是，高采樣1比特量化情況下的信號(hào)到達(dá)角高精度估計(jì)方法，具體包括以下步驟：

步驟1，產(chǎn)生由M個(gè)陣元構(gòu)成的均勻陣列，其僅接收的一個(gè)快拍信號(hào)y，該快拍信號(hào)經(jīng)過1比特量化采樣形成；

步驟2，建立以p范數(shù)為稀疏信號(hào)表達(dá)方式和外加抗噪聲干擾的不敏感支持向量機(jī)回歸的ε-SVR模型，利用該非凸優(yōu)化模型構(gòu)建原始信號(hào)稀疏表示方式；

步驟3，運(yùn)用基于交替方向乘子法(ADMM)求解該非凸優(yōu)化模型，求出稀疏表示系數(shù)，由稀疏表示系數(shù)確定信號(hào)的波達(dá)方向。

本發(fā)明的有益效果是，本發(fā)明高采樣1比特量化情況下的信號(hào)到達(dá)角高精度估計(jì)方法，能在采樣量少(如單快拍)，無關(guān)信源是否相干，不需計(jì)算自相關(guān)矩陣的情況下，高精度估計(jì)信號(hào)波達(dá)方向，并且基于1比特量化采樣情況下的信號(hào)波達(dá)方向估計(jì)更易于在工程中實(shí)現(xiàn)。

附圖說明

圖1是本發(fā)明用于接收單快拍信號(hào)的均勻陣列結(jié)構(gòu)圖；

圖2是本發(fā)明方法得到的實(shí)際信號(hào)波達(dá)方向的估計(jì)圖。

具體實(shí)施方式

下面結(jié)合附圖和具體實(shí)施方式對(duì)本發(fā)明進(jìn)行詳細(xì)說明。

本發(fā)明高采樣1比特量化情況下的信號(hào)到達(dá)角高精度估計(jì)方法，具體包括以下步驟：

步驟1，產(chǎn)生由M個(gè)陣元構(gòu)成的均勻陣列，其僅接收的一個(gè)快拍信號(hào)y，該快拍信號(hào)經(jīng)過1比特量化采樣形成；

假設(shè)接收信號(hào)滿足窄帶信號(hào)，即信號(hào)經(jīng)過陣列長(zhǎng)度的時(shí)間遠(yuǎn)小于信號(hào)的相干時(shí)間，如圖1所示，來波方向?yàn)棣?sub>i(i＝1,2,....I)，入射N根天線，陣元間隔為d的均勻線陣的陣列響應(yīng)矢量為：

定義方向矩陣為：

相對(duì)應(yīng)的信號(hào)幅度為則經(jīng)采樣1比特量化后的信號(hào)y：

其中V是高斯白噪聲,Tr是給定閾值，若測(cè)量值實(shí)部(虛部)大于該閾值實(shí)部(虛部)，則y的實(shí)部(虛部)為1，否則為-1。

步驟2，建立以p范數(shù)為稀疏信號(hào)表達(dá)方式和外加抗噪聲干擾的不敏感支持向量機(jī)回歸的ε-SVR模型，利用該非凸優(yōu)化模型構(gòu)建原始信號(hào)稀疏表示方式；

考慮這N個(gè)信號(hào)的入射角取值范圍為：[-90 90]，為此，將[-90 90]均勻細(xì)分為181個(gè)均勻空間并建立導(dǎo)向矢量過完備字典：

A＝[a(θ₁) a(θ₂) … a(θ₁₈₁)]∈C^M×181 (4)

相對(duì)應(yīng)的信號(hào)幅度S＝(s₁,s₂,…s₁₈₁)；為了從單快拍y中估計(jì)θ_i，建立以p范數(shù)為稀疏約束項(xiàng)和不敏感ε-SVR關(guān)于角度估計(jì)的稀疏表示模型：

(·)_n,r和(·)_n,i分別表示·的第n分量的實(shí)部和虛部，ζ_n,ζ^*_n,ε>0。

步驟3，運(yùn)用基于交替方向乘子法(ADMM)求解該非凸優(yōu)化模型，求出稀疏表示系數(shù)，由稀疏表示系數(shù)確定信號(hào)的波達(dá)方向。

(5)式是一個(gè)非凸優(yōu)化模型，使用交替方向乘子法(ADMM)迭代算法求解在過完備字典A(基)下的稀疏表示系數(shù)S，從而得出信號(hào)的波達(dá)方向。具體步驟如下：

首先，根據(jù)ADMM，引入輔助變量t＝S，構(gòu)造如下的拉格朗日函數(shù)：

關(guān)于該模型p,ρ,c,ε的選取很重要，懲罰因子c>0，控制對(duì)錯(cuò)分樣本的懲罰程度，c值過大過小都會(huì)影響系統(tǒng)的泛化能力；不敏感損失系數(shù)ε過小，則回歸精度高，但支持向量增多，反之，回歸精度降低；p＝1，則演變?yōu)镾的1范數(shù)，p＝0，則演變?yōu)镾的0范數(shù)，0范數(shù)表示向量非零元素的個(gè)數(shù)，實(shí)際上0范數(shù)不可解，所以，為了更稀疏約束信號(hào)，選取0<p<1，p越接近0，S越稀疏；ρ控制||S||^p和的均衡，以及控制ADMM算法迭代步長(zhǎng)；

然后，分布優(yōu)化如下：

①固定t,λ，求解S，則(6)式等價(jià)為：

為進(jìn)一步簡(jiǎn)化(7)，即：

其中，I＝181，顯然(8)式取最小值時(shí)，可以得出S與同相位，所以直接考慮復(fù)數(shù)模的運(yùn)算，并把(8)式轉(zhuǎn)化為I個(gè)子函數(shù)求解問題：

其中是t,S,λ的模，

因?yàn)樵谏戏謩e單調(diào)遞增和遞減，所以的全局最小值出現(xiàn)在上，即(9)式等價(jià)于：

現(xiàn)在討論在上的凹凸性，并分別對(duì)求1階導(dǎo)，2階導(dǎo)，3階導(dǎo)：

對(duì)p分三種情況討論(10)式的最小值點(diǎn)：

(1)1<p≤2

因?yàn)闀r(shí)，故是凸的，并且所以在區(qū)間有唯一解，使用二分法對(duì)(11)求零點(diǎn)即為的最小值點(diǎn)；

(2)p＝1

當(dāng)p＝1時(shí)(10)式等價(jià)于：

則的最小值點(diǎn)為：

(3)0<p<1

0<p<1時(shí)故是單增的，的凹凸性取決于的符號(hào)，記的解為：

當(dāng)是凹的，故的最小值點(diǎn)在0或之間；

當(dāng)可能是正的或負(fù)的，分和兩個(gè)區(qū)間討論：

i.在上，是凹的，故的最小值點(diǎn)在0或之間；

ii.在上，是凸的，且如果是單增的，故上的局部最小值點(diǎn)如果使用二分法求取上的局部最小值點(diǎn)；

顯然，的全局最小值點(diǎn)通過比較和上的局部最小值而獲得；

②固定S,λ，求取t，則(6)式等價(jià)為：

目標(biāo)函數(shù)進(jìn)一步化簡(jiǎn)為：

其中令w＝t-r，則(17)式等價(jià)于：

其中進(jìn)一步令real(·)和imag(·)分別表示·的實(shí)部與虛部，則(19)化簡(jiǎn)為：

其中，H＝[real(A)-imag(A)]∈C^M×362，H_n∈C^1×362為H的第n行，R＝[imag(A)real(A)]∈C^M×362，R_n∈C^1×362為R的第n行，D＝real(Ar+V-Tr)∈C^M×1，D1＝imag(Ar+V-Tr)∈C^M×1；

利用拉格朗日乘子法來解決這個(gè)約束優(yōu)化問題，構(gòu)造如下的拉格朗日函數(shù)：

其中，是拉格朗日乘子；

根據(jù)最優(yōu)化理論，L分別對(duì)W,ξ,ξ^*求偏微分，并令所得等式為0，得到下式：

將(22)式帶入(21)式得到對(duì)偶優(yōu)化問題：

通過矩陣簡(jiǎn)化，可以轉(zhuǎn)化(23)式為：

其中Q＝[H^T*diag(real(y)R^T*diag(imag(y)]∈C^362×2M，

P＝[D^T*diag(real(y)D1^T*diag(imag(y)]∈C^1×2M，T＝[real(y)^T imag(y)^T]

公式(24)用CVX優(yōu)化工具包，求解出β，進(jìn)而解出：

W＝Qβ (25)

由可以得出w，利用w＝t-r，求出：

t(k+1)＝w+r (26)

k為迭代次數(shù)；

③固定S,t，求λ，即：

λ(k+1)＝λ(k)+ρ(t-S) (27)

最后，①，②，③步循環(huán)迭代一定次數(shù)，直到S,λ,t收斂，即求得信號(hào)的波達(dá)方向的稀疏表示系數(shù)S，然后作出信號(hào)角度頻譜圖，即S的幅值曲線，尖峰處角度，即為估計(jì)的信號(hào)到達(dá)角。

實(shí)施例1

步驟1：有M＝20個(gè)陣元，其僅接受一個(gè)快拍信號(hào)y，噪聲方差為0.01，信噪比為16.99db，設(shè)定閾值Tr為接受信號(hào)的最大值與最小值之間的隨機(jī)值。假設(shè)信號(hào)元來波方向?yàn)?20°,-40°，經(jīng)過(3)式1比特采樣量化得到y(tǒng)。

步驟2：建立以p范數(shù)為稀疏約束項(xiàng)和不敏感ε-SVR關(guān)于角度估計(jì)的稀疏表示模型：

(·)_n,r和(·)_n,i分別表示·的第n分量的實(shí)部和虛部，

步驟3：運(yùn)用基于交替方向乘子法(ADMM)求解該非凸優(yōu)化模型，求出稀疏表示系數(shù)，由稀疏表示系數(shù)確定信號(hào)的波達(dá)方向。

根據(jù)ADMM，引入輔助變量t＝S,構(gòu)造如下的拉格朗日函數(shù)：

因?yàn)锳DMM算法是迭代算法，所以首先初始化所求參數(shù)t,λ，然后分布優(yōu)化如下：

①固定t,λ，求解S，則(29)式等價(jià)為：

為進(jìn)一步簡(jiǎn)化(30)，即：

其中，I＝181，顯然(31)式取最小值時(shí)，可以得出S與同相位，所以直接考慮復(fù)數(shù)模的運(yùn)算，并把(31)式轉(zhuǎn)化為I個(gè)子函數(shù)求解問題：

其中是t,S,λ的模，

因?yàn)樵谏戏謩e單調(diào)遞增和遞減，所以的全局最小值出現(xiàn)在上，即(32)式等價(jià)于：

現(xiàn)在討論在上的凹凸性，并分別對(duì)求1階導(dǎo)，2階導(dǎo)，3階導(dǎo)：

對(duì)p分三種情況討論(33)式的最小值點(diǎn)：

(1)1<p≤2

因?yàn)闀r(shí)，故是凸的，并且所以在區(qū)間有唯一解，使用二分法對(duì)(34)求零點(diǎn)即為的最小值點(diǎn)。

(2)p＝1

當(dāng)p＝1時(shí)(33)式等價(jià)于：

則的最小值點(diǎn)為：

(3)0<p<1

0<p<1時(shí)故是單增的，的凹凸性取決于的符號(hào)，記的解為：

當(dāng)是凹的，故的最小值點(diǎn)在0或之間。

當(dāng)可能是正的或負(fù)的，分和兩個(gè)區(qū)間討論：

i.在上，是凹的，故的最小值點(diǎn)在0或之間。

ii.在上,是凸的，且如果是單增的，故上的局部最小值點(diǎn)如果使用二分法求取上的局部最小值點(diǎn)。

顯然，的全局最小值點(diǎn)通過比較和上的局部最小值而獲得。

②固定S,λ，求取t,則(29)式等價(jià)為：

目標(biāo)函數(shù)進(jìn)一步化簡(jiǎn)為：

其中令w＝t-r，則(40)式等價(jià)于：

其中進(jìn)一步令real(·)和imag(·)分別表示·的實(shí)部與虛部，則(42)化簡(jiǎn)為：

其中，H＝[real(A)-imag(A)]∈C^20×362，H_n∈C^1×362為H的第n行，R＝[imag(A)real(A)]∈C^20×362，R_n∈C^1×362為R的第n行，D＝real(Ar+V-Tr)∈C^20×1，D1＝imag(Ar+V-Tr)∈C^20×1。

利用拉格朗日乘子法來解決這個(gè)約束優(yōu)化問題，構(gòu)造如下的拉格朗日函數(shù)：

其中，是拉格朗日乘子。

根據(jù)最優(yōu)化理論，L分別對(duì)W,ξ,ξ^*求偏微分，并令所得等式為0，得到下式：

將(45)式帶入(44)式得到對(duì)偶優(yōu)化問題：

通過矩陣簡(jiǎn)化，可以轉(zhuǎn)化(46)式為：

其中Q＝[H^T*diag(real(y)R^T*diag(imag(y)]∈C^362×40，

P＝[D^T*diag(real(y)D1^T*diag(imag(y)]∈C^1×40，T＝[real(y)^T imag(y)^T]

公式(47)用CVX優(yōu)化工具包，可以求解出β，進(jìn)而解出：

W＝Qβ (48)

由可以得出w，利用w＝t-r，求出：

t(k+1)＝w+r (49)

③固定S,t，求λ，即：

λ(k+1)＝λ(k)+ρ(t-S) (50)

①，②，③步循環(huán)迭代一定次數(shù)k，直到S,λ,t收斂，即求得S，

作出S曲線，尖峰點(diǎn)處的角度即所對(duì)應(yīng)的信號(hào)到達(dá)角。如圖2所示，信號(hào)到達(dá)角為-20°,-40°，高精度估計(jì)了這兩個(gè)信源方向。