專利名稱：一種用于線列陣的高精度測向方法
技術領域：
本發(fā)明涉及聲納信號處理技術領域，更具體地，本發(fā)明涉及一種用于線 列陣的高精度測向方法。
背景技術：
聲納信號處理的 一個基本問題就是測向，以確定感興趣的空間信號的方 向。高精度測向技術，是指采用某種算法，對陣列接受到的目標輻射或反射 的時空二維信號進行處理，從而得到目標方位的高精度估計。高精度測向技 術在聲納、雷達等領域有著廣泛的應用。
在目前典型的用于線列陣的高精度測向算法包括"J. Capon， High—Resolution Frequency—Wavenumber Spectrum Analysis, Proceedings of the IEEE, Vol. 57， No. 8, August 1969"中的Capon法，"孫長瑜等,二維 頻域波束形成方法，應用聲學，1995"中的窄帶二維頻域波束形成法，
"J"ames. A. Cadzow， Direction—of Arrival Estimation Using Signal Subspace Modeling, IEEE transaction on aerospace and electronic systems, Vol. 25, No. 1, January 1992"中的多重信號分類(MUSIC)法，"Shahram Shahbazpanahi, Robust Adaptive Beamforming for General-Rank Signal Models, IEEE transactions of signal processing，vol.51,no.9,September 2003，，中的自適應波束形成法以及"Roy R等，ESPRIT—A Subspace Rotation approach to estimation of parameters of cissoids in noise, IEEE Trans. ASSP, 1986"中的旋轉不變量信號參數估計方法(ESPRIT)法。
Capon法，是指通過保持目標信號位置輸出功率不變，最小化系統總輸出
5功率，在目標信號位置形成極大響應，其它方位輸出最小。表示為
^n(wH*R*W)，wH*a = 1,其中，a為指定方向的方向矢量，R是陣列輸出協方 差矩陣，w為加權向量。對于此方法，如果輸入信號的長度為N，則相應的要 計算N*N維矩陣的逆矩陣，隨著N的增長，計算量會急劇增加，對整個系統 的計算能力要求極高，而一般的線列陣的陣元個數N是幾十乃至幾百，如此 巨大的矩陣求逆導致該方法實用性差。
窄帶二維頻域波束形成法在窄帶目標信號輸入的條件下，在時間域和空 間域對輸入陣列數據進行快速傅立葉變換(FFT),在頻域得到波束形成。該 方法只能處理窄帶信號，與測向時經常遇到的寬帶信號不匹配，實用性差。
MUSIC法，是一種典型的空間語估計測向方法，首先估計出信號源的個數， 其次對信號輸出協方差矩陣進行特征值分解，最后對信號方位進行掃描，找 到對應信號源個數的譜峰。該方法對陣列模型失配十分敏感，存在陣列誤差 時其估計性能明顯下降。此方法要求信號源是非相關的，而且噪聲為具有相
同功率的空間白噪聲，這兩點在實際應用中一般難以滿足，導致系統性能不 穩(wěn)定，系統穩(wěn)健性差，另外該方法需要對信號輸出協方差矩陣進行特征值分 解，涉及到大矩陣求逆的問題，也帶來計算量巨大的問題。
自適應波束形成法，是指系統根據環(huán)境噪聲場的變化，不斷的自動調節(jié) 本身的參數以適應周圍環(huán)境，抑制干擾并檢測出有用信號。該方法的系統性 能穩(wěn)定性較差，而且不斷調整參數也帶來計算量大的問題，實際應用困難。
ESPRIT法，是指利用陣列傳感器分組所得的信號子空間中隱含的旋轉不 變量實現測量來波方向的方法。子空間的不變量通過建立自相關和互相關矩 陣來實現。將原陣列分解為兩個子陣，通過最小二乘(Least Square, LS)或 總體最小二乘(Total Least Square, TLS)擬合來估計子陣列的相移矩陣，得 到目標信號方位估計。在計算中涉及元素間互相關值的估計時，從自相關和 互相關矩陣中扣除估得的噪聲方差有時卻很致命，并使整體效果不佳，使得
6系統的穩(wěn)定性差，實用性能較差。
總體而言，由于水聲工作環(huán)境的惡劣性(例如，噪聲^^莫型易失配、海洋 傳播信道的復雜、陣形易失配等)，以上方法在現階段的實際應用中，或者 計算量巨大，系統難以負荷，或者實用性差，穩(wěn)健性不高，難以在實際的應 用中得到采用。

發(fā)明內容
為克服現有高精度側向技術中的實用性差、穩(wěn)健性不高并且難以處理寬 帶目標信號的缺陷，本發(fā)明提供一種用于線列陣的高精度測向方法。
根據本發(fā)明的一個方面，提出了一種用于線列陣的高精度測向方法，包
括
步驟IO)、對所述線列陣接收的時空二維信號，進行頻域傳統波束形成， 在波束域內得到粗測目標角度；
步驟20)、將所述線列陣劃分為若干個子陣，在所述粗測目標角度內， 對各子陣接收的時空二維信號分別進行頻域傳統波束形成并聚焦形成頻域-
波數數據；
步驟30 )、對所述頻域-波數數據在所述粗測目標角度范圍內根據寬帶聚 焦最小方差無失真方法進行處理，獲取目標信號方位。 其中，步驟10)進一步包括
步驟IIO)、獲取各陣元的所述時空二維信號的頻域數據，對所述頻域數 據進行空間域補零，并對每一個頻點在空間域上做快速傅立葉變換和i普移， 得到頻點和波束對應數據；
步驟120)、在所述每一個頻點，依據頻率-波數網格對所述頻點和波束 對應數據進行校正，得到頻率-波數數據矩陣；
步驟130)、對所述頻率-波數數據矩陣在頻域做反快速傅立葉變換，并 進行波束形成處理，獲得各搜索方位上的輸出功率；步驟140)、確定大于"&定閾值的各方位的輸出功率的信號為具有所述粗 測目標角度的信號。
其中，步驟20)包括
步驟210)、將線列陣在邏輯上劃分為多個子陣，使用所述各子陣接收時 空二維信號；
步驟220 )、對所述時空二維信號在時間域上做快速傅立葉變換，獲得各 陣元的頻域數據，并且生成在粗測角度范圍和處理帶寬內各個頻率分量上的 子陣間延時表；
步驟230 )、對子陣的頻域數據，在空間域上補零，并在目標頻段內的每 一個頻率分量上做空間域的快速傅立葉變換和譜移，得到頻點和波束對應數 據；
步驟240 )、在處理帶寬內的每一個頻點上，對所述頻點和波束對應^:據 進行校正，得到頻率-波數數據矩陣，通過反快速傅立葉變換，得到時間-波 數數據矩陣。
其中，步驟30)包括
步驟310)、在所述粗測目標角度范圍的每一個角度，將所述時間-波數 數據矩陣重組為處理帶寬x子陣數維的數據矩陣；
步驟320 )、計算所述處理帶寬x子陣數維的數據矩陣的協方差矩陣及其 逆矩陣，并且根據導向矢量為全'T，矢量的寬帶聚焦最小方差無失真算法， 獲得粗測角度范圍內的輸出功率；
步驟330 )、搜索所述輸出功率，將大于設定閾值的信號確定為目標信號， 輸出所述目標信號方位。
其中，步驟110)還包括對所述時空二維信號在時間域上做快速傅立葉 變換，獲得各陣元的頻域數據。
其中，步驟IIO)中，所述空間域補零可以在原凝::悟后直接補零，也可以 在原數據中插值補零。其中，步驟210)中，對線列陣在邏輯上進行劃分的方式包括子陣不重 疊、子陣重疊和子陣交叉。
其中，步驟230 )中，所述空間域補零可以在原^:據后直接補零，也可以 在原數據中插值補零。
其中，步驟240 )還包括在處理帶寬的每一個頻點上，依據頻率-波數 網格和所述子陣間延時表對所述二維信號校正，得到頻率-波數數據矩陣。
其中，對于所述接收到的時空二維信號的長度較長，分為多個長度較短 的信號。
通過應用本發(fā)明的方法，采用頻域波束形成，充分利用快速傅立葉變換， 計算速度快，便于DSP工程實現，能實現實時處理；采用快速傅立葉變換， 在頻域上運算，將寬帶信號分解為多個頻點信號，對每一個頻點進行處理， 適應于寬帶噪聲目標測向；采用計算效率高的頻域常規(guī)波束形成進行目標方 位的粗測，然后采用降低維數及計算復雜度的高精度分子陣的MVDR算法進行 目標方位精測，算法穩(wěn)健性高，目標探測精度高。


圖l是現有技術中的常規(guī)線列陣示意圖2是高精度測向方法的總體流程圖3是線列陣子陣劃分示意圖4是本發(fā)明高精度測向方法的具體流程圖5是頻率-波數網格的示意圖6是根據一個實施例的聲納數據處理和顯控系統。
具體實施例方式
下面結合附圖和具體實施例對本發(fā)明提供的一種用于線列陣的高精度測 向方法進4于詳細描述。圖1示出現有技術中的常規(guī)線列陣，如圖1所示，用于接收空間信號的
線列陣102由32個常規(guī)的無指向性水聽器101組成，水聽器101的接收信號 中心頻率為5000KHz，任意兩個水聽器之間間隔半個波長，線列陣012的聲振 段的總長度是4. 65米，線列陣102可以安裝在潛艇或無人潛航器上。
在根據本發(fā)明的實施例所提供的方法中，如圖2所示，通過利用快速傅 立葉變換計算的高效性，對整個線列陣接收到的時空二維數據在頻域實現快 速的常規(guī)波束形成，在波束域掃描后得到寬帶目標信號的粗測方位；然后在 粗測方位附近利用寬帶聚焦MVDR算法進行寬帶目標方位的高精度估計，獲取 高精度目標方位估計輸出?？紤]到常規(guī)MVDR算法涉及到計算N*N維矩陣的逆 矩陣的問題，將線列陣在邏輯上進行劃分為若干個子陣，在每一個子陣內進 行快速的頻域常規(guī)波束形成，將每一個子陣視為一個虛擬陣元，再進行寬帶 聚焦MVDR算法，大大降低了計算的復雜性，同時也能獲得高精度的寬帶目標 方位估計。
線列陣的子陣劃分方式有很多種，其中，有三種常見的子陣劃分方式， 如圖3所示，方式201為子陣不重疊的劃分方式，若各個子陣中心之間的距 離過大，會導致出現柵瓣，但對于不很長的線列陣，不會出現該問題；方式 202為子陣重疊的劃分方式，由于陣元的重疊，可以避免柵瓣的出現，但相應 增加了計算量；方式203為子陣交叉的劃分方式，雖然保留整條線陣的分辨 率，但是如果子陣數過多，將出現柵瓣現象，在以下根據本發(fā)明的實施例中， 采用201子陣不重疊的劃分方式。本領域內的普通技術人員可以理解，使用 其它兩種劃分方式也同樣可以實現本發(fā)明的方法，并達到所要求的技術效果。
在根據本發(fā)明的一個實施例中，圖4詳細示出一種用于線列陣的高精度 測向方法的流程。
如圖4所示，對單個線列陣接收到的時空二維信號，進行頻域傳統波束 形成(傳統波束形成延時相加波束形成，頻域傳統波束形成是將延時相加波 束形成方法在頻域實現)在波束域內掃描得到粗測的目標角度。步驟501:用單個線列陣接收時空二維信號，從該線列陣的各陣元獲得時
域數據。
步驟502:對線列陣中的每個陣元，分別耳又一段接收數據在時間域上做快 速傅立葉變換(FFT)，得到各陣元的頻域數據，組成一個頻域-空間域數據 矩陣。
本步驟中，需要注意，因為取的數據長度越長，即數據矩陣維數越大， 則運算速度會相應減慢。為了保證運算速度能夠滿足實時處理的要求，數據 長度不宜過大。 一般FFT運算點數在1024點及以下均能滿足要求。
步驟503:在空間域上對步驟502獲取的凄丈據補零，然后對處理帶寬內的 每一個頻點在空間域上做快速傅立葉變換，得到頻域-波數矩陣；將獲得數據 的零頻分量移至譜中心，得到一組與各頻點和波束相對應的數據，即頻率-波 數數據矩陣；
由于線列陣總陣元數M —般不大，為了在空間域上做FFT后頻率-波數網 格(如圖5所示)中有更"稠密"的分布，4吏得樣本點能接近任意的導向角，
在空間域上對數據補零至^點，^遠遠大于M。
空間域補零可以有多種選擇，既可以在原lt據后直接補零，也可以在原 數據中插值補零，只要補零后的數據長度滿足要求即可。
本領域的普通技術人員可以理解，將零頻分量移至譜中心是一種常用的 處理方法，在許多文獻中都有說明。
步驟504:為以下所估l的校正，需要預先產生頻率-波數網格。
具體地說，本步驟中，預先產生的頻率-波數網格由以下各式定義 "A,附)=row"d(iVp *厶*r(em) + 0.5) r(ej = c/*sin(6>m)/c.
， ， 5
其中，*是頻點，^，—表示頻率-波數網格中頻點^和波數附位置的數據，
厶人分別為采樣頻率和第k個頻點頻率，"是陣元間距，^是掃描角，^是 空間域補零后點數，W是FFT點數，陽朋"W表示取最接近x的整數，x(n,:)表示取矩陣x的第n行。
步驟505:在每一個頻點上依據預先產生的頻率-波數網格對由步驟503 所得的頻域-波數矩陣進行校正。
采用FFT實現頻域波束形成所對應的角度是波數域上的值，不是真正的 角度，為了得到對應真實的角度的值， 一般需要對其進行插值，插值一般意 味著巨大的計算量，本發(fā)明采用預先生成用來4交正的頻率-波數網格的方法進 行校正,該方法可參考文南史"Brian Maranda , Efficient digital beamforming in the frequency domain, 1989, L Acoustical Society of America"。
具體地說，本步驟中，頻率-波數網格才交正由下式定義
，，:)=+ Wp /2)
其中，*是頻點，尸是空間域補零后FFT得到的數據，《是經頻率—波數
網格校正后數據，^是空間域補零后點數，x(n，)表示取矩陣x的第n行。
步驟506:對步驟505獲得的頻率-波數凄t據矩陣在頻域做反快速傅立葉 變換，獲得時間-波數數據矩陣；
步驟507:將步驟506獲得的時間-波數數據矩陣的每一列的共軛轉置與 其自身相乘，得到在搜索方位上的輸出功率；搜索各方位輸出，大于設定閾 值的，認為是目標信號位置，輸出這些目標信號方位；本步驟是波束形成步 驟，時間-波數數據矩陣的一列表示一個搜索方位， 一行表示一個時間點的信
—;一
搜索角度的步長越小，輸出的波束越多，探測的精度也越高。輸出的目 標信號方位估計即為粗測的目標方位，這些^f直用來確定下一步寬帶聚焦MVDR 算法處理的角度范圍。
對單個線列陣在邏輯上進行分割，分為上面結合圖3所述的若干個子陣， 在粗測目標角度范圍內，對各子陣的時空二維信號分別進行頻域傳統波束形 成，由于維數下降，降低了計算量，隨后進行聚焦(延時各子陣波束形成結 果)送往下一步；步驟508:對單個線列陣在邏輯上進行分割，分割為N—sub個子陣，用各 子陣接收時空二維信號；所劃分的子陣參見圖3所示；
步驟509 ):根據步驟507 )獲得的粗測的目標方位，確定需要精測的角 度區(qū)域；
本步驟中，精測區(qū)域的確定，要根據實際的項目要求和硬件條件來確定。 本實施例中取粗測角度左右各3度為需要精測的角度區(qū)域，這也要求粗測的 目標方位的誤差不能超出± 3度， 一般而言頻域常M^波束形成能達到這個要求。
步驟510):根據步驟508 )確定的子陣劃分方案和步驟509 )確定的精 測角度范圍，生成子陣間的延時表；
具體的說，本步驟中子陣間延時表的生成由下式決定
fife—0,力=2tt * z * ^ * ~2-^
其中，,為采樣頻率，"是陣元間距，Z是處理帶寬內的頻點位置，7'是 精測區(qū)域內掃描角度，^是FFT的點數，M。是子陣內的陣元數。
步驟511):對所劃分的每一個子陣進行頻域常規(guī)波束形成；
本步驟中，對子陣數據進行頻域常規(guī)波束形成與步驟502 )至步驟506 ) 基本一致，區(qū)別在于在進行步驟505 )的處理時，進行校正時需要同時乘以一 個延時量，該延時量由子陣序號與子陣延時表的乘積決定；另外各子陣數據 的處理后的時間-波數數據矩陣均要存儲；
具體的說，本步驟中，頻率-波數網格校正后數據乘以一個延時量由下式 定義
=尸(&,"&,:)+ /2)， 尸3:)=尸2:) * exp() ** ^to;;(A:,:))
其中，^是子陣序號，*是頻點，尸是空間域補零后FFT得到的數據，A是經頻率-波數網格校正后數據，《是經頻率-波數網格校正后數據乘以一個
延時量后的數據，^是空間域補零后點數，x (n，)表示取矩陣x的第n行。
步驟512):在粗測目標角度范圍內，在每一個角度上，對上述步驟輸出 的各子陣的時間-波數數據矩陣，重組為頻率-子陣數據矩陣；計算獲得的頻 率-子陣數據矩陣的協方差矩陣，得到N-sub*N_sub的協方差矩陣；N—sub為 子陣數；計算該協方差矩陣的逆矩陣；
513) 按照導向矢量為全'T，矢量的寬帶聚焦MVDR算法公式，在粗測角 度范圍內的每一個角度上，計算系統功率輸出。
本步驟中，全1導向矢量的出現是因為在前述的子陣頻域常規(guī)波束形成 中已經考慮了各陣元間、各子陣間的延時，所以此時的導向矢量無需在考慮 延時問題，故取為全1的導向矢量，即意味著聲波到達各陣元間是沒有延時 差的。
具體的說，本步驟中的導向矢量為全"1"矢量的寬帶聚焦MVDR算法公 式為
s = [1,1,1,......l]， N—sub個l;
pMVDR=abs( *丄r)
s承R 承s ,
其中R"為步驟512)輸出的協方差矩陣的逆矩陣；PMVDR為輸出功率。
514) 搜索步驟513)的輸出，大于設定閾值的，即為目標信號位置，輸出 這些精測的目標信號方位。
對長度很大的數據，可以采取分段處理的方法，對每個分段數據進行上 述用于線列陣的快速穩(wěn)健的高精度測向，再積累輸出各個分段的檢測結果。
圖6示出一種可以應用本發(fā)明所述方法的現有的聲納數據處理和顯控系 統。其中，拖曳線列陣102和前置電路601從每個水聽器接收到模擬信號。 前置電路601包括前放、濾波和其它常規(guī)的電^各。每個通道的模擬信號輸入 到A/D轉換器602得到數字信號。從A/D轉換器602出來的是多通道的數字
14數據流，每個通道的數據流對應于一個水聽器接收到的模擬信號。將這些數
據流輸入到微型處理器603。經微型處理器603處理后的輸出信息可以存儲在 數據存儲器件605，比如磁盤存儲設備中，或直接輸出到顯示設備606上顯示。
微型處理器603首先將接收到的數據流存儲到動態(tài)存取區(qū)604,在輸入滿 足處理要求數量的數據流后就開始處理。執(zhí)行上述圖4的流程，包括時間域 FFT 5 02,空間域補零后做FFT,并移零頻分量至鐠中心503,頻率-波數校正 505，頻率-波數網格504,頻域IFFT 506，波束形成507，產生子陣間延時表 510，子陣頻域常規(guī)波束形成511,子陣聚焦，求相關矩陣及其逆矩陣512， 全'T，導向矢量的寬帶聚焦MVDR算法513及目標方位搜索514。這些程序儲 存在動態(tài)存取區(qū)604中。
因為本發(fā)明中的A/D轉換器602輸出的是多通道的數據流，因此可以采 用多片的微型處理器來并行處理。能夠實現圖6功能的其它一些硬件設備， 比如專用硬件、基于應用的集成電路(ASIC) 、 DSP等都可以用來代替微型處 理器603。
最后應說明的是，以上實施例僅用以描述本發(fā)明的技術方案而不是對本 技術方法進行限制，本發(fā)明在應用上可以延伸為其他的修改、變化、應用和 實施例，并且因此認為所有這樣的修改、變化、應用、實施例都在本發(fā)明的 精神和教導范圍內。
1權利要求
1、一種用于線列陣的高精度測向方法，包括步驟10)、對所述線列陣接收的時空二維信號，進行頻域傳統波束形成，在波束域內得到粗測目標角度；步驟20)、將所述線列陣劃分為若干個子陣，在所述粗測目標角度內，對各子陣接收的時空二維信號分別進行頻域傳統波束形成并聚焦形成頻域-波數數據；步驟30)、對所述頻域-波數數據在所述粗測目標角度范圍內根據寬帶聚焦最小方差無失真方法進行處理，獲取目標信號方位。
2、 權利要求l的方法，其中，步驟10)進一步包括步驟IIO)、獲取各陣元的所述時空二維信號的頻域數據，對所述頻域數 據進行空間域補零，并對每一個頻點在空間域上做快速傅立葉變換和譜移， 得到頻點和波束對應lt據；步驟120)、在所述每一個頻點，依據頻率-波數網檔_對所述頻點和波束 對應數據進行校正，得到頻率-波數數據矩陣；步驟130)、對所述頻率-波數數據矩陣在頻域做反快速傅立葉變換，并 進行波束形成處理，獲得各搜索方位上的輸出功率；步驟140)、確定大于設定閾值的各方位的輸出功率的信號為具有所述粗 測目標角度的信號。
3、 權利要求2的方法，其中，步驟110)還包括對所述時空二維信號 在時間域上做快速傅立葉變換，獲得各陣元的頻域數據。
4、 權利要求2的方法，其中，步驟IIO)中，所述空間域補零可以在原 數據后直接補零，也可以在原數據中插值補零。
5、 權利要求l的方法，其中，步驟20)進一步包括步驟21Q)、將線列陣在邏輯上劃分為多個子陣，使用所述各子陣接收時空二維信號；步驟220 )、對所述時空二維信號在時間域上做快速傅立葉變換，獲得各陣元的頻域數據，并且生成在粗測角度范圍和處理帶寬內各個頻率分量上的子陣間延時表；步驟230 )、對子陣的頻域數據，在空間域上補零，并在目標頻段內的每 一個頻率分量上做空間域的快速傅立葉變換和譜移，得到頻點和波束對應數 據；步驟240 )、在處理帶寬內的每一個頻點上，對所述頻點和波束對應數據 進行校正，得到頻率-波數數據矩陣，通過反快速傅立葉變換，得到時間-波 數數據矩陣。
6、 權利要求5的方法，其中，步驟210)中，對線列陣在邏輯上進行劃 分的方式包括子陣不重疊、子陣重疊和子陣交叉。
7、 權利要求5的方法，其中，步驟230 )中，所述空間域補零可以在原 數據后直接補零，也可以在原數據中插值補零。
8、 權利要求5的方法，其中，步驟240 )還包括在處理帶寬內的每一 個頻點上，依據頻率-波數網格和所述子陣間延時表對所述二維信號進行校 正，得到頻率-波數數據矩陣。
9、 權利要求5的方法，其中，步驟30)進一步包括步驟310)、在所述粗測目標角度范圍的每一個角度，將所述時間-波數 數據矩陣重組為處理帶寬x子陣數維的數據矩陣；步驟320 )、計算所述處理帶寬x子陣數維的數據矩陣的協方差矩陣及其 逆矩陣，并且才艮據導向矢量為全'T，矢量的寬帶聚焦最小方差無失真算法， 獲得粗測角度范圍內的輸出功率；步驟330 )、搜索所述輸出功率，將大于設定閾值的信號確定為目標信號， 輸出所述目標信號方位。
10、 權利要求1的方法，其中，對于長度較長的所述接收到的時空二維信號，將其分為長度較短的多個信號。
全文摘要
本發(fā)明公開一種用于線列陣的高精度測向方法，包括對線列陣接收的時空二維信號，進行頻域傳統波束形成，在波束域內掃描得到粗測目標角度；將所述線列陣劃分為若干個子陣，在所述粗測目標角度范圍內，對各子陣接收的時空二維信號分別進行頻域傳統波束形成并聚焦形成頻域-波數數據；對所述頻域-波數數據在所述粗測目標角度范圍內根據寬帶聚焦最小方差無失真方法進行處理，獲取目標信號方位。通過應用本發(fā)明，計算速度快，能實現實時處理，算法穩(wěn)健性高，目標探測精度高。
文檔編號G01S3/80GK101470187SQ20081014681
公開日2009年7月1日 申請日期2008年8月25日 優(yōu)先權日2007年12月26日
發(fā)明者葉青華, 宇 李, 彪 田, 黃海寧 申請人:中國科學院聲學研究所