0 引言

    信源定位是陣列信號處理的一個重要研究方向，現代科學技術正在蓬勃發展，定位技術更是廣泛應用在航空航天、交通、勘探、導航等領域^[1-2]。例如，軍工應用中，雷達需要定位技術作為新的補充，只需利用傳感器接收到的信號即可確定目標的位置，而且在定位過程中受到的干擾和攻擊可降低到最低水平。因此，信源定位技術因其諸多優點及廣泛的應用前景成為學者們的研究重點。

    現在常用的信源定位方法主要分成三類：(1)利用時延估計進行定位^[3]，先進行陣元間時間差的計算，再根據陣列結構與信源模型估計出信源的位置。該定位方法計算量較小，易于實時實現，但僅在單信源定位中得到了廣泛的應用。(2)利用波束形成法進行定位^[4]，無需計算時間差，需要通過優化目標函數來實現信源定位，但是實際應用環境中會出現多個最優值即偽峰，所以峰值搜索過程較為復雜。(3)利用高分辨空間譜估計進行定位^[5-6]，可同時高分辨率地定位多個聲源，但是對實驗環境要求較高。目前大多定位算法都是以窄帶信號為研究對象提出的。寬帶信號中含有較大帶寬和多個中頻，直接使用窄帶信號對應算法時造成的相位差無法忽略，因此針對寬帶信號中的信號頻率也有處理方法。文獻[7-8]已將最大似然法（ML）擴展至寬帶信號處理中，但是由于涉及到高度非線性化及非凸性和最大化似然函數，其計算復雜度較高。另一類是Kaveh在1985年提出的具有高分辨DOA估計的相干信號子空間法(CSSM)[9-11]，通過引入“聚焦”思想，使得不同頻率上的觀測值全部變換到確定頻率上，后對各聚焦后信號的協方差矩陣進行平均，得到聚焦后的協方差矩陣，最后獲得寬帶信號的高分辨DOA估計，但該方法估計性能易受信源方位預估精度的影響。本文提出了一種針對寬帶信源定位的聯合質心收縮與聯合可控響應功率和相位變換(SRP-PHAT)的信源定位方法^[4，12]，根據對寬帶信源獲得的DOA估計，利用本文提出的算法構建出信源所在的初始區域，通過質心收縮、布置虛擬信源點等方法逐步收縮區域，最終獲得信源定位估計。

1 分布式陣列模型與信號模型

    本文提出的寬帶信源定位方法采用的分布式陣列系統由四個均勻線陣組成，如圖1所示。

    K個均勻線陣任意分布在分布式系統中，放置于水平面上，且每個線陣含有N個傳感器，第i個信源所在位置坐標為(x_si，y_si，z_si)，第k個線陣中心坐標為(x_k，y_k，z_k)，第k個線陣對第i個信源的方位角為θ_ik，由于各線陣是任意分布，故線陣中心連線圍成一任意多邊形區域V，信源可分布在該多邊形區域內的任意位置。

導向矢量α_k(f_j，θ_M)可表示為：

    假設噪聲為高斯白噪聲，且與各信號均相互獨立。陣列流型矩陣可表示為：

2 信源定位方法

2.1 DOA估計

    DOA估計是信源定位方法的關鍵步驟，是信源定位確定初始區域的關鍵。傳統的MUSIC算法對常見的陣列都普遍適用，故傳統MUSIC算法的空間譜估計公式可定義為：

式中U_N則是N×(N-M)噪聲子空間。

    由于寬帶信號中有多個中心頻率，故需將不同中心頻率的信號子空間映射到一個最佳的聚焦頻率f₀，通過式（8）選取最佳聚焦頻率：

    通過譜峰搜索在空間譜P_MUSIC(f₀，θ)中得到。但是，對于密集信源來說，使用傳統MUSIC算法進行DOA估計性能較差，而且在少量傳感器時甚至無法區分出信源的數目。為提高密集信源的DOA估計精度，本文引入了MUSIC空間譜的群延遲函數來提高密集聲源在少量傳感器時DOA估計的精度。Group Delay函數利用了MUSIC空間譜相位信息的負差分形式，用其與MUSIC空間譜函數以乘積的形式來消除群延遲譜的偽峰，提高DOA估計精度。Group Delay譜定義為：

其中Φ(f₀，θ)為MUSIC空間譜相位信息，表示為：

    對P_MGD(f₀，θ)進行譜峰搜索后可得到精確的DOA估計值。

2.2 構建初始區域及質心收縮

    根據式(11)得到的各子陣的DOA估計值、第k個子陣中心的坐標A_k=(x_k，y_k)和第i個信源對第k個子陣的DOA估計值θ_ki，利用點斜式方程，做出一條子陣中心出發的平面射線。所有子陣做出的直線圍成一個初始區域，各子陣中心連線與各射線交于一點B_ki=(x_ki，y_ki)，則各射線的交點I_kp=(I_x，I_y)即A_k B_ki與A_p B_pi的交點可通過式(12)、式(13)計算為：

    所有的交點I_kp均保存在集合I中，所有的交點的橫縱坐標的最大值和最小值可構成初始區域V₁。通過I中所有交點的坐標可以計算出兩兩交點之間的相對距離，并找出最大相對距離D_max，而且D_max是隨集合I的變化而變化。

    根據I中的交點，計算初始區域的質心：

其中|I|代表集合I中交點的數目。同時可計算出各交點到區域質心的歐幾里得距離，表示為：

其中s代表維數，D={D_k}。根據得到的歐幾里得距離，找出最小距離對應交點將其插入到新集合P_j中，并從集合I中刪除。然后再計算集合P_j中交點所圍成區域的質心Cent(P_j)，同時計算出集合I中剩余交點到質心Cent(P_j)的歐幾里得距離：

根據上述得到的最大相對距離D_max和集合I、P_j，同時找出集合D′中的最小值min(D′)。將min(D′)與D_max進行比較，當min(D′)≤D_max時，將D′中最小距離對應交點插入到集合P_j中，同時將該點從集合I中刪除。重復上述質心收縮步驟，直到min(D′)>D_max時，停止迭代。最后，更新后的集合P_j中的所有交點的橫縱坐標最值[x_max  x_min  y_max  y_min]圍成粗略收縮區域V₂。

2.3 最大互功率譜收縮定位法

    本節介紹利用SRP-PHAT在區域V2中隨機布置虛擬信源點進行區域收縮。首先將聯合可控功率函數定義為：

τ(x，k)為信號從信源沿著得到的DOA估計值到第k個線陣的到達時間長度。由于各虛擬信源點的坐標均已知，可通過虛擬信源點到第k個線陣中心的距離d_k和信號傳播速度c獲得時間延遲τ(x，k)，即：

    本文的信源定位算法由DOA估計、區域質心收縮和最大互功率譜收縮組成。定位算法步驟如下：

    (1)由式（11）獲得的DOA估計值和各線陣中心得到各線陣到信源的射線。

    (2)根據式（12）、式（13）獲得各射線直接的交點。

    (3)根據式（14）、式（15）和min(D′)≤Dmax獲得最終更新后的交點集合P_j。

    (4)P_j中交點圍成了區域V₂。

    (5)在V₂中隨機布置虛擬信源點，再根據式（16）計算各虛擬信源點到各線陣的互功率譜。

    (6)找出前E個最大互功率譜值對應點，創建新的收縮區域V₃。

    (7)重復步驟(5)和步驟(6)，直到滿足然后計算V_peak中虛擬信源點坐標平均值即得到信源定位估計值。

3 仿真實驗和分析

    實驗仿真證明了本文所提算法的良好性能，圖1為實驗采用的分布式陣列。各子陣均是采用四個陣元的均勻線陣，線陣中心連線圍成了一個4.6 m×5 m的矩形區域，且密集信源在該矩形區域內部，信源1坐標為(2.95 m，2.95 m)，信源2坐標為(1.2 m，3.45 m)，用*代表信源。圖中子陣1的中心坐標為A₁=(0 m，0 m)，子陣2、子陣3及子陣4的中心坐標分別為A₂=(5 m，0 m)、A₃=(5 m，4.6 m)、A₄=(0 m，4.6 m)。其中子陣1水平放置，并以其中心坐標為坐標原點，即α=0°；子陣2相對于子陣1逆時針旋轉了45°，即β=45°；子陣3相對于子陣1逆時針旋轉了90°，即γ=90°；子陣4相對于子陣1逆時針旋轉了30°，即η=30°。本文使用實際錄音信號為實驗的聲源信號,聲音持續時間為3 s。

    實驗對虛擬信源點數F取值為100，最大互功率譜值對應點數E取值為10。門限值δ取值為1%。根據密集信源位置，圖1中兩信源對應子陣1中心的角度值為35°和40°，利用CSSM和Group Delay譜的方法對這兩個信源進行DOA估計，其譜估計如圖3所示，所用噪聲為高斯白噪聲，信噪比SNR為5 dB。由圖2可知，本文使用的DOA估計算法分辨率較高。

    根據前面得到的密集信源的DOA估計值，以聲源1為例，利用各子陣中心點和對應DOA估計值作出相應射線，射線所圍成的區域即為初始搜索區域，如圖3所示。

    根據確定的初始區域，利用質心收縮和最大互功率譜收縮逐步縮小區域，圖4（a）是兩密集信源的區域收縮過程，其中實線代表信源1的收縮過程，虛線代表信源2的收縮過程。圖4（b）則是對圖4（a）中圈中區域的細節放大圖。圖中兩實心點為信源位置估計。

    為了對本文算法性能進行評估，提出了位置估計的均方差（RMSE）隨信噪比變化的評估方法。本次實驗共做200次蒙特卡羅實驗，其中RMSE定義為：

4 結束語

    本文著重討論了二維空間中寬帶信號的信源高精度DOA估計和信源定位方法。首先，利用Group Delay函數和CSSM算法對傳統MUSIC算法的空間譜加以優化，消除無用偽峰，提高寬帶信號的DOA估計精度。利用獲得的DOA估計值確定初始搜索區域，然后在二維平面中使用質心收縮和最大互功率譜收縮法對初始區域逐步收縮，得到最終信源的定位估計。本文提出的算法計算簡單便捷，僅需DOA一個參數即可獲得信源的定位估計。仿真實驗結果證明了本文算法的DOA估計精度較高，SNR對信源定位估計的RMSE的影響較小。但是本文算法較為依賴DOA估計確定的初始區域，DOA估計精度的大小直接影響到算法的運算量，可以考慮引入其他參數來控制DOA精度導致的算法計算量過大的問題。

參考文獻

[1] 殷晴青，何培宇.相干語音信號的兩陣元DOA估計新算法[J].信號處理，2013，29(9)：1169-1175.

[2] 寧強，方向，潘俊，等.圓形傳感器陣列多運動聲源二維波達方向估計[J].解放軍理工大學學報：自然科學版，2013，14(3)：283-287.

[3] OMOLOGO M，SVAIZER P.Use of the crosspower-spectrum phase in acoustic event location[J].IEEE Transactions on Speech & Audio Processing，1997，5(3)：288-292.

[4] DIBIASE J H.A high-accuracy, low-latency technique for talker localization in reverberant environments[D].Brown University，2000.

[5] 居太亮，彭啟琮，邵懷宗，等.基于任意麥克風陣列的聲源二維DOA估計算法研究[J].通信學報，2005，26(8)：129-133.

[6] CHEN J C，YAO K，HUDSON R E.Acoustic source localization and beamforming：theory and practice[J].Eurasip Journal on Advances in Signal Processing，2003(4)：359-370.

[7] DORON M，WEISS A J，MESSER H.Maximum-likelihood direction finding of wide-band sources[J].IEEE Transactions on Signal Processing，1993，41(1)：411-414.

[8] CHEN J C，HUDSON R E，YAO K.Maximum-likelihood source localization and unknown sensor location estimation for wideband signals in the near-field[J].IEEE Transactions on Signal Processing，2002，50(8)：1843-1854.

[9] KHAN M A，RAO R K，WANG X.Performance of quadratic and exponential multiuser chirp spread spectrum communication systems[C].Performance Evaluation of Computer and Telecommunication Systems(SPECTS)，2013 International Symposium on.IEEE，2013：58-63.

[10] YOON Y S，KAPLAN L M，MCCLELLAN J H.TOPS：new DOA estimator for wideband signals[J].IEEE Transactions on Signal Processing，2006，54(6)：1977-1989.

[11] 張小飛，汪飛，徐大專.陣列信號處理的理論和應用[M].長沙：國防工業出版社，2010.

[12] ASTAPOV S，BERDNIKOVA J，PREDEN J S.A method of initial search region reduction for acoustic localization in distributed systems[C].Mixed Design of Integrated Circuits and Systems(MIXDES)，2013 Proceedings of the 20th International Conference.IEEE，2013：451-456.

[13] KUMAR L，TRIPATHI A，HEGDE R M.Robust multisource localization over planar arrays using music-group delay spectrum[J].IEEE Transactions on Signal Processing，2014，62(17)：4627-4636.