基于稀疏矩陣填充的級聯毫米波雷達高分辨測角方法

摘 要： 級聯毫米波雷達將多個射頻芯片進行級聯，可實現大規模多輸入多輸出（multiple-input multiple-output， MIMO）天線陣列。然而，現有毫米波雷達的測角方法存在精確度不足、旁瓣過高等問題。對此，提出一種基于稀疏矩陣填充的級聯毫米波雷達高分辨測角方法。首先，建立級聯毫米波雷達二維稀疏陣列的陣列信號模型。然后，將該陣列的單快拍信號構造為低秩漢克爾（Hankel）矩陣，使用加權非精確增廣拉格朗日乘子法進行稀疏矩陣填充，并將填充后的矩陣轉化為均勻虛擬陣列的單快拍信號。最后，使用常規波束形成（conventional beamforming， CBF）方法估計目標角度信息并進行配對。仿真和實測實驗結果表明，所提方法能夠精確地估計目標角度信息，并有效抑制峰值旁瓣電平。

關鍵詞： 級聯毫米波雷達; 波達方向估計; 稀疏矩陣; 矩陣填充

中圖分類號： TN 958

文獻標志碼： A

DOI：10.12305/j.issn.1001-506X.2024.08.11

High resolution angle measurement method for cascaded millimeter wave

radar based on sparse matrix completion

ZHANG Minglong¹， WU Yulin¹， WEI Wenqiang¹， SHEN Yuanjie¹，

GUO Shisheng^1，2， CUI Guolong^1，2，*

（1. School of Information and Communication Engineering， University of Electronic Science and

Technology of China， Chengdu 611731， China; 2. Yangtze Delta Region Institute （Quzhou），

University of Electronic Science and Technology of China， Quzhou 324003， China）

Abstract： Cascaded millimeter wave radar cascades multiple radio frequency chips to realize large-scale multiple-input multiple-output （MIMO） antenna array. However， the existing angle measurement methods of millimeter wave radar have problems such as insufficient accuracy and high sidelobes. To solve these problems， a high resolution angle measurement method for cascaded millimeter wave radar is proposed based on sparse matrix completion. Firstly， a two-dimensional sparse array signal model of the cascaded millimeter wave radar is established. Then， the single snapshot signal of the array is converted into a low-rank Hankel matrix， the weighted inexact augmented Lagrange multiplier method is used for sparse matrix completion， and the completed matrix is converted into the single snapshot signal of uniform virtual array. Finally， conventional beamforming （CBF） method is used to estimate the angles of the targets and the angles of the targets are paired. The simulation and experimental results show that the proposed method can accurately estimate the angles of the targets and effectively suppress peak sidelobe level.

Keywords： cascaded millimeter wave radar; direction of arrival （DOA） estimation; sparse matrix; matrix completion

0 引 言

隨著自動駕駛技術的發展，毫米波雷達已經成為高級駕駛輔助系統（advanced driving assistance system， ADAS）不可或缺的一部分^［1］。毫米波雷達是一種利用毫米波頻段（30～300 GHz）的電磁波進行探測的雷達。相比于其他ADAS傳感器，如激光雷達、攝像頭等，其具有全天時、全天候的探測能力，可在惡劣天氣和低光照條件下進行有效的目標探測和跟蹤^［2-3］。同時，毫米波雷達還具有硬件成本低、體積小等特點，因此毫米波雷達在自動駕駛領域得到了廣泛應用，如自適應巡航控制、緊急制動、泊車輔助等^［4-5］。

在毫米波雷達信號處理過程中，波達方向（direction of arrival， DOA）估計是重要步驟之一，影響著毫米波雷達對目標的定位精度^［6］。現有的毫米波雷達DOA估計方法包括波束形成類^［7-9］、子空間類^［10-11］、迭代自適應方法（iterative adaptive approach， IAA）類^［12-15］、壓縮感知類^［16］、深度學習類^［17］方法等。文獻［18］基于單片微波集成電路（mono-lithic microwave integrated circuit， MMIC）芯片毫米波雷達，利用數字波束形成方法估計目標的方位角和仰角，但是該方法的角度估計精度有限。文獻［19］基于單芯片毫米波雷達，通過最小方差無失真響應（minimum variance distortionless response， MVDR）方法進行目標測角，以消除多徑假目標;文獻［20］使用多重信號分類（multiple signal classification， MUSIC）方法估計目標角度，通過單芯片毫米波雷達實現目標檢測，但是上述兩種方法的角度估計性能依賴于多快拍的陣列協方差矩陣。文獻［21］基于單芯片毫米波雷達，提出了一種三維卷積神經網絡以用于DOA估計，不過該方法需要大量的數據進行訓練，較難應用于實際工程中。但是，單芯片毫米波雷達收發陣元數量較少，角度分辨率有限，并且其性能容易受到信號衰減、多徑效應等影響，而級聯毫米波雷達將多個MMIC芯片進行級聯，可實現大規模多輸入多輸出（multiple-input multiple-output， MIMO）天線陣列，提高角度分辨率、抗干擾能力和探測精度^［22-24］。對于級聯毫米波雷達均勻陣列，文獻［25］使用MUSIC方法分析級聯毫米波雷達的角度估計性能，通過對比實驗證明其角度估計精度優于單芯片毫米波雷達。文獻［26］結合先驗信息和IAA，提出了一種快速局部迭代搜索算法，使用級聯毫米波雷達在單快拍條件下估計目標角度。壓縮感知方法在文獻［27］中被用于級聯毫米波雷達目標角度估計，以實現對周圍環境的準確感知。雖然上述的IAA和壓縮感知方法能夠基于單快拍信號得出角度估計結果，但是需要較高的計算成本和信噪比條件，并且網格失配情況會降低壓縮感知方法角度估計的精度。

級聯毫米波雷達的角度分辨率取決于天線陣列的虛擬孔徑。在確定收發天線陣元數量的前提下，使用傳統均勻布陣方法得到的均勻陣列虛擬孔徑較小、角度分辨率較低。而采用稀疏布陣方法得到的稀疏陣列能夠在不增加天線陣元和硬件成本的條件下，獲得更大的陣列虛擬孔徑，實現更高的角度分辨率，從而提高級聯毫米波雷達的目標檢測性能，在有限的天線陣元資源下，挖掘更多的潛能^［28-30］。然而，上述DOA估計方法在被應用于稀疏陣列時存在精度不足、高旁瓣等問題。在級聯毫米波雷達稀疏陣列的基礎上，可通過矩陣填充（matrix completion， MC）方法恢復接收信號矩陣中的缺失數據，再進行DOA估計。該方法可提高缺失陣元數據時的DOA估計精度，同時避免網格失配問題^［31-34］，為實現高精度雷達成像和三維環境感知提供支撐。傳統的基于核范數最小化的MC算法包括奇異值閾值算法^［35］、加速近端梯度算法^［36］、交替方向乘子法^［37］、增廣拉格朗日乘子法^［38］等方法。為提高核范數最小化方法的靈活性，文獻［39］指定合理的核范數權重，實現加權核范數最小化，提高估計損壞的潛在輸入數據的能力;文獻［40］提出截斷核范數方法，通過截斷較大奇異值之和，更好地解決低秩矩陣逼近問題。文獻［41］使用奇異值閾值算法恢復陣列信號協方差矩陣中的缺失元素，利用旋轉不變技術估計目標角度參數。文獻［42］結合加速近端梯度算法和MUSIC算法以恢復接收信號矩陣和估計DOA^［42］。文獻［43］基于交替方向乘子法和聯合重加權以填充稀疏快拍矩陣，通過快速傅里葉變換（fast Fourier transform， FFT）獲得目標方位角和仰角。文獻［44］將一維稀疏陣列的快拍信號構造為漢克爾（Hankel）矩陣，使用奇異值閾值算法填充信號數據矩陣，再通過FFT算法實現一維DOA估計;文獻［45］將數據排列成具有多重低秩Hankel結構的低秩增強形式，通過核范數最小化方法實現MC，恢復缺失數據;文獻［46］為減輕MIMO雷達天線陣元故障帶來的負面影響，將四重Hankel矩陣的結構化處理與MC方法相結合，有效地恢復陣列信號協方差矩陣的行和列的連續損壞數據，提高天線陣列的角度估計性能^［46］。然而，上述方法只通過理論分析與數字仿真驗證了算法的有效性，其算法性能在實際場景中有待驗證。本文將MC引入到級聯毫米波雷達DOA估計，借助Hankel矩陣的特殊構型，通過迭代的方式恢復稀疏虛擬陣列信號以提高DOA估計性能，并通過數字仿真與實驗數據驗證了算法的有效性。

本文提出了一種基于稀疏MC的級聯毫米波雷達高分辨測角方法，與已有研究相比，本文的主要貢獻包括：① 構建了級聯毫米波雷達二維MIMO天線陣列的陣列回波信號模型;② 將陣列的單快拍信號構造為低秩Hankel矩陣，使用加權非精確增廣拉格朗日乘子法進行稀疏MC，獲取均勻虛擬陣列的單快拍信號，通過常規波束形成（conventional beamforming， CBF）方法估計目標方位角和仰角信息，并進行配對;③ 通過仿真實驗與實際場景測試，驗證了所提方法的有效性。結果表明，所提方法能夠準確估計目標的角度信息，并有效抑制旁瓣對估計結果的影響，可幫助級聯毫米波雷達實現多目標識別、點云成像等任務。

1 信號模型

假設級聯毫米波雷達由M個MMIC射頻芯片級聯組成，每個毫米波雷達芯片具有N_t個發射陣元和N_r個接收陣元。系統采用時分復用（time division multiplexing， TDM）技術，發射陣元發射調頻連續波（frequency modulated continuous wave， FMCW）信號。所有發射陣元和接收陣元構成級聯毫米波雷達二維MIMO天線陣列，在三維空間中對目標方位角、俯仰角進行估計。

建立如圖1所示的級聯毫米波雷達陣列模型，其中級聯天線陣列位于x-o-z平面，紅色表示發射陣列，由MN_t個發射陣元組成，藍色表示接收陣列，由MN_r個接收陣元組成。

發射陣元和接收陣元的x軸坐標可分別表示為

t_x=［t_x，1，t_x，2，…，t_x，MNt］^T（1）

r_x=［r_x，1，r_x，2，…，r_x，MNr］^T（2）

發射陣元和接收陣元的y軸坐標均為0，z軸坐標可分別表示為

t_z=［t_z，1，t_z，2，…，t_z，MNt］^T（3）

r_z=［r_z，1，r_z，2，…，r_z，MNr］^T（4）

假設遠場一目標位于坐標（x_o，y_o，z_o）處，方位角與仰角分別為φ和θ，以徑向速度V運動，則目標與第i個發射陣元、第j個接收陣元之間的距離分別為

R_t，i=（x_o－t_x（i））²+y²_o+（z_o－t_z（i））²（5）

R_r，j=（x_o－r_x（j））²+y²_o+（z_o－r_z（j））²（6）

發射陣元發射的Chirp信號為

s（t）=exp（j2πf₀t+jπKt²）， t∈［0，T］（7）

式中：f₀為載波頻率;K=B/T為調頻斜率，B為信號帶寬;T為信號時寬。對于級聯毫米波雷達第i個發射陣元所發射的第n個Chirp信號，n=1，2，…，N，N表示發射的Chirp信號數目，第j個接收陣元所接收到的回波信號可表示為

s_i，j（t，n）=s（t－（τ_i，j+2V（t+ΔT_n）/c））+v_i，j（t，n）=

exp（j2πf₀（t－（τ_i，j+2V（t+ΔT_n）/c）））·

s^b_i，j（t，n）+v_i，j（t，n）（8）

s^b_i，j（t，n）=exp（jπK（t－（τ_i，j+2V（t+ΔT_n）/c））²）（9）

式中：v_i，j（t，n）表示均值為0、方差為σ²的復高斯白噪聲;ΔT_n=（n－1）T;c表示光速;電磁傳播時延τ_i，j=（R_t，i+R_r，j）/c?？紤]到2V（t+ΔT_n）/c與τ_i，j相比是一個非常小的值，故忽略s^b_i，j（t，n）中包含的2V（t+ΔT_n）/c項，接收信號與發射信號經過混頻后，中頻信號可表示為

x_i，j（t，n）=exp（j2πf₀τ_i，j）x^b_i，j（t，n）+v_i，j（t，n） （10）

x^b_i，j（t，n）=exp（j2π（Kτ_i，j+f_d）t－jπKτ²_i，j+j2πf_dΔT）（11）

式中：f_d=2V/λ表示多普勒頻率;λ為信號波長。由于系統采用TDM技術，當目標運動時，不同發射陣元形成的虛擬陣列之間存在多普勒相位偏移。相鄰發射順序的發射陣元之間的多普勒相位偏移為

Δ?_d=4πVTλ（12）

則中頻信號為

x_i，j（t，n）=exp（j（2πf₀τ_i，j+（i－1）Δ?_d））·

x^b_i，j（t，n）+v_i，j（t，n）（13）

考慮到目標位于遠場，將式（5）和式（6）分別重寫為

R_t，i≈R_t，1+（t_x（i）－t_x（1））sin φcos θ+

（t_z（i）－t_z（1））sin θ（14）

R_r，j≈R_r，1+（r_x（j）－r_x（1））sin φcos θ+

（r_z（j）－r_z（1））sin θ（15）

故目標與第i個發射陣元、第j個接收陣元之間的電磁傳播時延可表示為

τ_i，j=R_t，i+R_r，jc=τ_1，1+Δτ_i，j（16）

式中：Δτ_i，j為第i個發射陣元、第j個接收陣元與第1個發射陣元、第1個接收陣元之間的波程差。由于其是非常小的值，故忽略x^b_i，j（t，n）中的Δτ_i，j項，將信號重寫為

x_i，j（t，n）=exp（j（2πf₀Δτ_i，j+（i－1）Δφ_d））·

x⁰（t，n）+v_i，j（t，n）（17）

式中：x⁰（t，n）=exp（j2πf₀τ_1，1）x^b_1，1（t，n）。

結合所有虛擬陣列的通道，將x_1，1（t，n），x_1，2（t，n），…，x_{MNt，MNr}（t，n）排列為向量x（t，n），則虛擬陣列接收信號為

x（t，n）=a（φ，θ）x⁰（t，n）+v（t，n）（18）

a（φ，θ），v（t，n）可分別表示為

a（φ，θ）=a_t（φ，θ）a_r（φ，θ） （19）

v（t，n）=［v_1，1（t，n），v_1，2（t，n），…，v_MNt，MNr（t，n）］^T（20）

式中：a（φ，θ）表示導向矢量;a_t（φ，θ）和a_r（φ，θ）分別為發射導向矢量和接收導向矢量;表示克羅內克積。a_t（φ，θ）和a_r（φ，θ）分別表示為

a_t（φ，θ）=［exp{j2π（t_x，1sin φcos θ+t_z，1sin θ）/λ}，

exp{j2π（t_x，2sin φcos θ+t_z，2sin θ）/λ+jΔ?_d}，…，

exp{j2π（t_x，MNtsin φcos θ+t_z，MNtsin θ）/λ+j（MN_t－1）Δ?_d}］^T（21）

a_r（φ，θ）=［exp{j2π（r_x，1sin φcos θ+r_z，1sin θ）/λ}，

exp{j2π（r_x，2sin φcos θ+r_z，2sin θ）/λ}，…，

exp{j2π（r_x，MNrsin φcos θ+r_z，MNrsin θ）λ}］^T（22）

經過采樣之后，得到離散接收信號x（h，n），可表示為

x（h，n）=a（φ，θ）x⁰（h，n）+v（h，n），h=1，2，…，H（23）

式中：H表示采樣點數。對每個虛擬陣列通道的信號進行H-點離散傅里葉變換（discrete Fourier transform， DFT），結果可表示為

x^r（h-，n）=a（φ，θ）x^r，0（h-，n）+v^r（h-，n）， h-=1，2，…，H-（24）

式中：x^r，0（h-，n）為x⁰（h，n）的DFT結果;v^r（h-，n）是噪聲的DFT結果。不難發現，每一個Chirp信號中，所有虛擬陣列通道的距離維DFT的目標峰值索引均為h_p。取x^r（h-，n）中索引為h_p的值，表示為

x^r（h_p，n）=a（φ，θ）x^r，0（h_p，n）+v^r（h_p，n）（25）

2 基于MC的高分辨測角方法

本節將介紹基于稀疏MC的級聯毫米波雷達高分辨測角方法。所提方法首先將MIMO天線陣列的單快拍信號構造為低秩Hankel矩陣，然后使用加權非精確增廣拉格朗日乘子法進行稀疏MC，并將填充后的矩陣轉化為均勻虛擬陣列的單快拍信號，最后使用CBF方法估計目標的方位、俯仰角信息并進行配對。測角方法的具體流程圖如圖2所示。

2.1 稀疏MC

對于同一距離/多普勒單元內的少量目標，由于其空間散射點是稀疏的，故由虛擬陣列快拍信號構造的矩陣是低秩的。在結合稀疏與低秩先驗的基礎上，可將虛擬陣列的快拍信號構造為待填充的輸入矩陣^［31］。具體地，對所有虛擬陣列通道在目標距離單元的單快拍信號進行選取，以方位維度為例，將方位維度上選取的信號構造為低秩Hankel矩陣，作為待填充的矩陣。

虛擬陣列在目標距離單元的第n個快拍信號y（n）=x^r（h_p，n），其在方位維度選取后的信號為y_d（n）=［y_d（n，0），y_d（n，i_d），…，y_d（n，L）］^T，y_d（n，i_d）表示虛擬通道坐標為i_d的第n個快拍信號，i_d∈I_d，I_d為虛擬陣列通道坐標集合，相鄰陣元最小間距d=λ/2，Ld為虛擬孔徑。假設孔徑為Ld的均勻陣列在目標距離單元的第n個快拍信號為y_u（n）=［y_u（n，0），y_u（n，j_u），…，y_u（n，L）］^T，均勻陣列通道坐標j_u=0，1，2，…，L，y_u（n，j_u）可表示為

將均勻陣列快拍信號y_u（n）按照Hankel矩陣形式進行構造，得到的Hankel矩陣可以看作為無缺失元素的均勻陣列快拍信號的空間采樣矩陣，該矩陣具有低秩的先驗特性。通過求解Hankel矩陣核范數優化問題，可以較好地恢復均勻陣列快拍信號中的缺失元素^［47-48］。將該矩陣作為稀疏MC問題的輸入矩陣D，使用加權非精確增廣拉格朗日乘子法進行稀疏MC^［37-39］，填充后的目標矩陣A可求解如下：

minA|A|_w，*， A+E=D; π_Ω（E）=0（27）

式中：A_w，*=∑^r_l=1w_lσ_l（A）;σ_l（A）為矩陣A的第l個奇異值;w=［w₁，w₂，…，w_r］^T;w_l為σ_l（A）的非負權值;r為奇異值個數;E是A和D之間的誤差矩陣;Ω為D中非零項的索引集合;π_Ω是保持Ω中元素不變的線性運算符。

式（27）中，MC問題的部分增廣拉格朗日函數可表示為

L（A，E，Y，μ）=A_w，*+〈Y，D－A－E〉+μ2D－A－E2F（28）

式中：〈A，B〉表示矩陣A^TB的跡函數;Y為拉格朗日乘子;μgt;0為懲罰因子;·_F表示Frobenius范數。通過交替更新變量A，E，Y，μ以最小化部分增廣拉格朗日函數L（A，E，Y，μ）。具體來說，該問題可通過迭代以下步驟解決。

步驟 1 固定E_k、Y_k和μ_k不變，通過最小化L（A，E_k，Y_k，μ_k）以更新A_k+1，計算公式為

A_k+1=argminA L（A，E_k，Y_k，μ_k）=

argminAA_w，*+〈Y_k，D－A－E_k〉+μ_k2D－A－E_k2F=

argminAA_w，*+μ_k2D－A－E_k+1μ_kY_k2F（29）

迭代計算A_k+1時，通過指定合理的權重w，提高核范數最小化問題的靈活性，更好地恢復輸入矩陣中待填充的部分。由于輸入矩陣的較大奇異值量化了其基本主方向的信息，所以較大奇異值的收縮量較少，較小奇異值的收縮量較大^［39］。矩陣A進行奇異值分解后可表示為A=UΣ-V^T，Σ-=diag（［σ₁（A），σ₂（A），…，σ_r（A）］）。令Z=D-E_k+Y_k/μ_k，矩陣Z進行奇異值分解后的奇異值為Σ=diag（［σ₁（Z），σ₂（Z），…，σ_r（Z）］），σ_l（A）為

σ_l（A）=0， c₂lt;0

c₁+c₂2， c₂≥0（30）

式中：c₁=σ_l（Z）－ε;c₂=（σ_l（Z）+ε）²－4C;ε為極小的正值;C為正則化參數。

步驟 2 固定A_k+1、Y_k和μ_k不變，通過最小化L（A_k+1，E，Y_k，μ_k）以更新E_k+1，計算公式為

E_k+1=argminπΩ（E）=0L（A_k+1，E，Y_k，μ_k）=argminπΩ（E）=0A_k+1w，*+

〈Y_k，D－A_k+1－E〉+μ_k2D－A_k+1－E2F=

argminπΩ（E）=0μ_k2D－A_k+1－E+1μ_kY_k2F=

πD－A_k+1+1μ_kY_k（31）

更新矩陣E_k+1時，約束π_Ω（E）=0在最小化L（A_k+1，E，Y_k，μ_k）時強制執行。

步驟 3 固定A_k+1、E_k+1和μ_k不變，更新Y_k+1：

Y_k+1=Y_k+μ_k（D－A_k+1－E_k+1）（32）

步驟 4 更新懲罰因子μ_k+1：

μ_k+1=ρμ_k（33）

式中：ρ為更新因子。

迭代上述步驟直至滿足迭代終止條件或者達到最大迭代次數時，終止迭代運算，得到目標矩陣A。基于仿真、實測實驗的數據處理經驗，將最大迭代次數k_m設置為200，迭代終止條件設置為D－A_k－E_kF/D_Flt;ε₁，其中ε₁表示迭代終止閾值，ε₁=10^－3。取目標矩陣A的反對角線元素的平均值^［47-48］，可得到均勻虛擬陣列的單快拍信號y_o（n）。稀疏MC算法如算法1所示。

算法 1 稀疏CM算法輸入 輸入信號y_u（n），虛擬孔徑Ld，迭代次數k，最大迭代次數k_m

輸出 輸出信號y_o（n）

步驟 1 將y_u（n）按照Hankel矩陣形式進行構造，得到輸入矩陣D和索引集合Ω。

步驟 2

步驟 2.1 更新A

A_k+1=argminAL（A，E_k，Y_k，μ_k）

步驟 2.2 更新E

E_k+1=argminπΩ（E）=0L（A_k+1，E，Y_k，μ_k）

步驟 2.3 更新Y

Y_k+1=Y_k+μ_k（D－A_k+1－E_k+1）步驟 2.4 更新μ

μ_k+1=ρμ_k

步驟 2.5 k′=k+1

步驟 3 若未達到最大迭代次數k_m或迭代終止條件D－A_k－E_kF/D_Flt;ε₁，重復步驟2。

步驟 4 取目標矩陣A的反對角線元素的平均值，得到輸出信號y_o（n）。

2.2 角度估計

對于均勻虛擬陣列的單快拍信號y_o（n），通過CBF方法^［7-8］進行DOA估計。CBF空間譜為

P（φ）=b_w（φ）^HRb_w（φ）（34）

式中：R=y_o（n）y_o（n）^H;b_w（φ）由方位維度的導向矢量加窗后得到，本文使用漢明窗。通過稀疏CM獲取均勻虛擬陣列的單快拍信號，然后使用CBF方法估計空間譜，獲取目標的角度，本文將上述方法稱為矩陣填充常規波束形成（matrix completion conventional beam forming， MC-CBF）方法。相似地，在俯仰維度進行DOA估計之后，獲得目標俯角。之后將獲取的方位角、俯角信息進行配對，對于S個方位角不同的目標（s=1，2，…，S），構造代價函數：

V_s=argp=1，2，…，Pmaxa^H（φ_s，θ_p）E_nE^H_na（φ_s，θ_p）^－1（35）

式中：E_n表示由虛擬陣列通道的快拍信號所構造的噪聲子空間;φ_s和θ_p分別為獲取的目標方位角和俯角;P為獲取的俯角個數;·表示2范數。通過計算代價函數的最大值，可將目標與獲取的仰角信息進行配對，得到最終的角度估計結果。

3 實驗分析

本節通過級聯毫米波雷達仿真實驗和實際場景測試，驗證本文所提的MC-CBF方法的性能。將MC-CBF方法與幾種毫米波雷達DOA估計方法進行對比，包括CBF方法^［7-8］、MVDR方法^［9］、MUSIC方法^［10］、IAA方法^［12-14］，驗證本文所提方法的有效性和優越性。

3.1 仿真實驗

考慮一級聯毫米波雷達，設置為兩片四發四收MMIC級聯模式。發射陣元和接收陣元的排布如圖3所示，其中相鄰陣元間距均為半波長的整數倍，最小間距為半波長。發射陣元的x軸和z軸坐標分別為t_x=［0，8，24，28，32，9，10，11］^T、t_z=［0，0，0，0，0，1，4，6］^T，接收陣元的x軸和z軸坐標分別為r_x=［0，1，2，3，12，14，15，53］^T、r_z=［0，0，0，0，0，0，0，0］^T。系統采用TDM技術，發射陣元的發射順序id為e_t=［1，2，3，4，5，6，7，8］。

發射陣元與接收陣元形成的二維虛擬陣列示意圖如圖4所示。方位、俯仰維度的虛擬子陣列分別為Subarray₁和Subarray₂。方位稀疏虛擬子陣列Subarray₁具有40個通道，構成的虛擬孔徑為42.5λ，俯仰稀疏虛擬子陣列Subarray₂具有4個通道，構成的虛擬孔徑為3λ。

假設兩個目標距離均為5 m，速度均為0 m/s，目標1的（方位角，俯角）為（0°，0°），目標2的（方位角，俯角）為（10°，0°），接收回波信噪比為SNR=20 dB，載波頻率為77 GHz，信號帶寬為2.527 GHz，信號時寬為32 μs，采樣率為8 MHz，采樣點數為256，快拍數為1。

用不同方法對Subarray₁在目標距離單元的單快拍信號進行DOA估計，得到的空間譜對比結果如圖5所示。由圖5可知，對于Subarray₁的單快拍信號，所有方法的空間譜在目標方位角處均出現了譜峰。由于Subarray₁為稀疏陣列，直接使用CBF方法進行DOA估計時會出現較高的旁瓣電平。MVDR方法和MUSIC方法受到快拍數限制，無法準確構造陣列信號協方差矩陣，因此旁瓣很高。MC-CBF方法通過稀疏MC得到均勻虛擬陣列的單快拍信號后，再根據式（34）進行DOA估計，相比于CBF、MVDR、MUSIC、IAA方法，旁瓣電平得到了有效抑制。圖6為保持其他條件不變，將陣列接收回波的信噪比設置為SNR=－10 dB時，不同方法的空間譜對比結果。從圖6可以看出，MC-CBF方法在低信噪比條件下仍然具有更好的旁瓣抑制水平。

圖7為目標的距離方位角度譜，圖7（a）、圖7（b）、圖7（c）分別為對Subarray₁的單快拍信號使用CBF、IAA、MC-CBF方法得到的距離方位角度譜。從圖7可知，在原有稀疏虛擬陣列基礎上直接進行DOA估計，由于旁瓣影響，不易直接在方位維度檢測到兩個目標，而MC-CBF方法的距離方位角度譜具有兩個清晰的譜峰，對應著正確的距離和方位角度。在俯仰維度進行DOA估計之后，將獲取的方位角、俯角信息進行配對。圖8為三維空間定位結果，對比了目標在空間中的真實位置與使用MC-CBF方法的估計結果。觀察圖8可知，MC-CBF方法能夠準確地估計目標的角度信息。

通過蒙特卡羅實驗，采用均方根誤差（root mean square error， RMSE）衡量不同方法對于虛擬陣列的DOA估計精度，RMSE定義為

RMSE=1S1Q∑Ss=1∑Qq=1（θ-_s，q－θ_s）²（36）

式中：S為信源個數;Q為蒙特卡羅實驗次數。圖9為在蒙特卡羅實驗中，不同方法在不同信噪比條件下對方位維度進行DOA估計的RMSE曲線。將接收回波信噪比設置為－10 dB到20 dB之間，將步進值設置為2 dB。由圖9可見，所有方法的RMSE都隨著信噪比的增加而逐漸減小，而且MC-CBF方法具有更高的角度估計精度。

為進一步評價所提方法的旁瓣抑制性能，將接收回波的信噪比設置為－10 dB到20 dB之間，設置步進值為2 dB。圖10為不同信噪比條件下，不同方法對方位維度DOA估計的空間譜的峰值旁瓣電平（peak sidelobe level， PSL）對比結果。從圖10可以發現，在不同信噪比下，本文提出的方法均具有更好的旁瓣抑制效果。

3.2 實測實驗

為驗證所提方法的有效性以及對比不同方法在實際測量中的性能，使用德州儀器（Texas Instruments， TI）公司的四芯片AWR2243級聯毫米波雷達進行實測實驗。其發射與接收陣列如圖11所示，其中相鄰陣元間距均為半波長的整數倍，最小陣元間距為半波長。該級聯毫米波雷達由4個MMIC芯片級聯構成，共包含12個發射陣元和16個接收陣元，構成的二維虛擬陣列在方位、俯仰維度的虛擬孔徑分別為42.5λ和3λ。

該級聯毫米波雷達的發射陣元與接收陣元能夠覆蓋仿真部分中的兩片四發四收MMIC級聯毫米波雷達的發射陣元與接收陣元。在實測實驗中，通過合理使能對應的8個發射陣元和8個接收陣元，即可實現仿真部分的級聯毫米波雷達，形成相同的虛擬陣列。故TI級聯毫米波雷達使能發射陣元的x軸和z軸坐標分別為t_x=［0，8，24，28，32，9，10，11］^T、t_z=［0，0，0，0，0，1，4，6］^T，使能接收陣元的x軸和z軸坐標分別為r_x=［0，1，2，3，12，14，15，53］^T、r_z=［0，0，0，0，0，0，0，0］^T。同時，設置發射陣元的發射順序id為e_t=［12，10，6，5，4，3，2，1］，使發射陣元、接收陣元的x、z軸坐標，發射陣元的發射順序均與仿真部分的級聯毫米波雷達相同。

圖12為實測實驗場景，令兩個目標距離均為5 m，速度均為0 m/s，目標1的（方位角，俯角）為（－4°，－2°），目標2的（方位角，俯角）為（6°，－2°），載波頻率為77 GHz，信號帶寬為2.527 GHz，信號時寬為32 μs，采樣率為8 MHz，采樣點數為256，快拍數為1。

圖13為不同方法對方位維度的稀疏虛擬子陣列的單快拍信號進行DOA估計得到的空間譜。觀察圖13可知，對于Subarray₁的單快拍信號，所有方法的空間譜在目標方位角處均出現了譜峰，但是MC-CBF方法能夠獲得比其他方法更低的旁瓣電平。

圖14為實測實驗中的目標距離方位角度譜，圖14（a）、圖14（b）、圖14（c）分別為對方位維度的稀疏虛擬子陣列使用CBF方法、IAA方法、MC-CBF方法所得到的距離方位角度譜?？梢园l現，使用MC-CBF方法對方位維度稀疏虛擬子陣列的單快拍信號填充之后，再進行DOA估計，獲得的距離角度譜具有兩個更加清晰的譜峰，該方法抑制了旁瓣電平，能夠獲取精確的目標角度信息。在俯仰維度進行DOA估計之后，將目標的方位角和仰角進行配對，得到如圖15所示的三維空間實測定位結果。觀察圖15可知，在實測實驗中，MC-CBF方法能夠較為準確地估計目標的方位角、俯角信息。

為了避免隨機性的影響，對級聯毫米波雷達虛擬陣列接收信號的任意快拍信號進行處理，分析任意單快拍下，各方法對方位維度DOA估計的空間譜的旁瓣抑制性能。從圖16可以發現，在接收回波信噪比較高的條件下，對于任意單快拍，所提方法具有更好的旁瓣抑制效果。

4 結 論

本文提出了一種基于稀疏MC的級聯毫米波雷達高分辨測角方法。首先，建立了級聯毫米波雷達二維稀疏陣列的陣列信號模型。其次，將虛擬陣列的單快拍信號構造為低秩Hankel矩陣，將其作為待填充的矩陣。然后，基于該矩陣奇異值的稀疏特性，構建稀疏MC的優化模型，并使用加權非精確增廣拉格朗日乘子法求解待填充的矩陣。最后，將填充后的矩陣構造為均勻陣列的單快拍信號，使用CBF方法估計目標的角度信息并進行配對。仿真實驗和實測實驗證明了所提方法的可行性，結果表明該方法能夠精確地估計目標角度信息，而且可以有效抑制PSL。

參考文獻

［1］SUN S Q， PETROPULU A P， POOR H V. MIMO radar for advanced driver-assistance systems and autonomous driving： advantages and challenges［J］. IEEE Signal Processing Magazine， 2020， 37（4）： 98-117.

［2］HAKOBYAN G， YANG B. High-performance automotive radar： a review of signal processing algorithms and modulation schemes［J］. IEEE Signal Processing Magazine， 2019， 36（5）： 32-44.

［3］洪偉， 余超， 陳繼新， 等. 毫米波與太赫茲技術［J］. 中國科學： 信息科學， 2016， 46（8）： 1086-1107.

HONG W， YU C， CHEN J X， et al. Millimeter wave and terahertz technology［J］. Scientia Sinica Informationis， 2016， 46（8）： 1086-1107.

［4］PATOLE S M， TORLAK M， WANG D， et al. Automotive radars： a review of signal processing techniques［J］. IEEE Signal Processing Magazine， 2017， 34（2）： 22-35.

［5］WALDSCHMIDT C， HASCH J， MENZEL W. Automotive radar-from first efforts to future systems［J］. IEEE Journal of Microwaves， 2021， 1（1）： 135-148.

［6］WANG Z G， LI W H， PENG H L， et al. The development of compact patch antenna array with high angular resolution for 77 GHz FMCW radar applications［C］∥Proc.of the IEEE Electrical Design of Advanced Packaging and Systems， 2021.

［7］SZALAY Z， NAGY L. Target modeling， antenna array design and conventional beamforming algorithms for radar target DOA estimation［C］∥Proc.of the 17th International Conference on Transparent Optical Networks， 2015.

［8］王平. 汽車毫米波雷達目標DOA估計研究［D］. 成都： 電子科技大學， 2021.

WANG P. Research on target doa estimation of automotive millimeter wave radar［D］. Chengdu： University of Electronic Science and Technology of China， 2021.

［9］CAPON J. High-resolution frequency-wavenumber spectrum analysis［J］. Proceedings of the IEEE， 1969， 57（8）： 1408-1418.

［10］SCHMIDT R. Multiple emitter location and signal parameter estimation［J］. IEEE Trans.on Antennas and Propagation， 1986， 34（3）： 276-280.

［11］DONG Y Y， DONG C X， LIU W， et al. 2-D DOA estimation for L-shaped array with array aperture and snapshots extension techniques［J］. IEEE Signal Processing Letters， 2017， 24（4）： 495-499.

［12］ROBERTS W， STOICA P， LI J， et al. Iterative adaptive approaches to MIMO radar imaging［J］. IEEE Journal of Selected Topics in Signal Processing， 2010， 4（1）： 5-20.

［13］FORTUNATI S， GRASSO R， GINI F， et al. Single-snapshot DOA estimation using compressed sensing［C］∥Proc.of the IEEE International Conference on Acoustics， Speech and Signal Processing， 2014.

［14］趙佳楠. 車載毫米波雷達目標方位估計與識別技術研究［D］. 哈爾濱： 哈爾濱工業大學， 2021.

ZHAO J N. Research on target azimuth estimation and recognition technology of automotive millimeter wave radar［D］. Harbin： Harbin Institute of Technology， 2021.

［15］SUN W Z， SO H C， CHEN Y， et al. Approximate subspace-based iterative adaptive approach for fast two-dimensional spectral estimation［J］. IEEE Trans.on Signal Processing， 2014， 62（12）： 3220-3231.

［16］ROOS F， HUGLER P， TORRES L L T， et al. Compressed sensing based single snapshot DOA estimation for sparse MIMO radar arrays［C］∥Proc.of the 12th German Microwave Conference， 2019： 75-78.

［17］ANNALURU R S， MAZHER K U， HEATH R W. Deep learning based range and DOA estimation using low resolution FMCW radars［C］∥Proc.of the IEEE Statistical Signal Processing Workshop， 2021： 366-370.

［18］CHENG Y W， SU J R， CHEN H Y， et al. A new automotive radar 4D point clouds detector by using deep learning［C］∥Proc.of the IEEE International Conference on Acoustics， Speech and Signal Processing， 2021： 8398-8402.

［19］鄭晶月， 吳佩侖， 陳家輝， 等. 車載毫米波雷達多徑假目標分析與消除方法［J］. 系統工程與電子技術， 2024， 46（1）： 88-96.

ZHENG J Y， WU P L， CHEN J H， et al. Multipath 1 target analysis and elimination for automotive millimeter wave radar［J］. Systems Engineering and Electronics， 2024， 46（1）： 88-96.

［20］GAO X Y， XING G B， ROY S， et al. Experiments with mmwave automotive radar test-bed［C］∥Proc.of the 53rd Asilomar Conference on Signals， Systems， and Computers， 2019.

［21］SANG T H， CHIEN F T， CHANG C C， et al. DOA estimation for FMCW radar by 3D-CNN［J］. Sensors， 2021， 21（16）： 5319.

［22］GAO X Y， ROY S， XING G B， et al. Perception through 2D-MIMO FMCW automotive radar under adverse weather［C］∥Proc.of the IEEE International Conference on Autonomous Systems， 2021.

［23］YANIK M E， WANG D， TORLAK M. 3-D MIMO-SAR imaging using multi-chip cascaded millimeter-wave sensors［C］∥Proc.of the IEEE Global Conference on Signal and Information Processing， 2019.

［24］KITAMURA T， SUWA K. Doppler division multiplexed multiple-input multiple-output imaging using cascaded millimeter-wave radars［J］. IEEE Trans.on Microwave Theory and Techniques， 2022， 70（3）： 1571-1581.

［25］ZHAI X， JIN G， CAI Z R， et al. Performance analysis of DOA estimation for texas instruments mmwave radar sensors［C］∥Proc.of the Global Reliability and Prognostics and Health Mana-gement， 2022： 1-6.

［26］WEI W Q， LIU R T， YU X X， et al. Fast single-snapshot DOA estimation of coherent sources for distributed mmwave radar system［J］. IEEE Tran.on Circuits and Systems II： Express Briefs， 2022， 69（8）： 3615-3619.

［27］CORREAS-SERRANO A， GONZALEZ-HUICI M A. Experimental evaluation of compressive sensing for DOA estimation in automotive radar［C］∥Proc.of the 19th International Radar Symposium， 2018.

［28］SYEDA R Z， SAVELYEV T G， BEURDEN M C， et al. Sparse MIMO array for improved 3D mm-wave imaging radar［C］∥Proc.of the 17th European Radar Conference， 2021： 342-345.

［29］MA Y， MIAO C， ZHAO Y Y， et al. An MIMO radar system based on the sparse-array and its frequency migration calibration method［J］. Sensors， 2019， 19（16）： 3580.

［30］YANG F W， XU F， FIORANELLI F， et al. Practical investigation of a MIMO radar system capabilities for small drones detection［J］. IET Radar， Sonar amp; Navigation， 2021， 15（7）： 760-774.

［31］SUN S Q， BAJWA W U， PETROPULU A P. MIMO-MC radar： a MIMO radar approach based on matrix completion［J］. IEEE Trans.on Aerospace and Electronic Systems， 2015， 51（3）： 1839-1852.

［32］宋虎， 蔣迺倜， 劉溶， 等. 基于稀疏采樣陣列優化的APG-MUSIC算法［J］. 電子與信息學報， 2018， 40（6）： 1390-1396.

SONG H， JIANG N T， LIU L， et al. APG-MUSIC algorithm based on sparse sampling array optimization［J］. Journal of Electronics amp; Information Technology， 2018， 40（6）： 1390-1396.

［33］LIU S H， MAO Z H， ZHANG Y M D， et al. Rank minimization-based toeplitz reconstruction for DOA estimation using coprime array［J］. IEEE Communications Letters， 2021， 25（7）： 2265-2269.

［34］孫兵， 阮懷林， 吳晨曦， 等. 幅度相位誤差條件下的互質陣列DOA估計方法［J］. 系統工程與電子技術， 2021， 43（12）： 3488-3494.

SUN B， RUAN H L， WU C X， et al. DOA estimation method for coprime array under gain and phase error［J］. Systems Engineering and Electronics， 2021， 43（12）： 3488-3494.

［35］CAI J F， CANDES E J， SHEN Z W. A singular value thres-holding algorithm for matrix completion［J］. SIAM Journal on Optimization， 2010， 20（4）： 1956-1982.

［36］TOH K C， YUN S. An accelerated proximal gradient algorithm for nuclear norm regularized linear least squares problems［J］. Pacific Journal of Optimization， 2010， 6（3）： 615-640.

［37］BOYD S， PARIKH N， CHU E， et al. Distributed optimization and statistical learning via the alternating direction method of multipliers［J］. Foundations and Trends in Machine Learning， 2011， 3（1）： 1-122.

［38］ROCKAFELLAR R T. Augmented lagrange multiplier functions and duality in nonconvex programming［J］. SIAM Journal on Control， 1974， 12（2）： 268-285.

［39］GU S H， XIE Q， MENG D Y， et al. Weighted nuclear norm minimization and its applications to low level vision［J］. International Journal of Computer Vision， 2017， 121（2）： 183-208.

［40］HU Y， ZHANG D B， YE J P， et al. Fast and accurate matrix completion via truncated nuclear norm regularization［J］. IEEE Trans.on Pattern Analysis and Machine Intelligence， 2013， 35（9）： 2117-2130.

［41］XU H Y， ZHANG Y K， BA B， et al. Two-dimensional direction-of-arrival fast estimation of multiple signals with matrix completion theory in coprime planar array［J］. Sensors， 2018， 18（6）： 1741.

［42］LI M X， HU W， DI J Y， et al. An APG-MUSIC algorithm based-on optimized sampling array［C］∥Proc.of the MATEC Web of Conferences， 2018.

［43］DING J R， WANG Z Y， MA W D， et al. TDM-MIMO automotive radar point-cloud detection based on the 2-D hybrid sparse antenna array［J］. IEEE Trans.on Geoscience and Remote Sensing， 2022， 60： 5112215.

［44］SUN S Q， PETROPULU A P. A sparse linear array approach in automotive radars using matrix completion［C］∥Proc.of the IEEE International Conference on Acoustics， Speech and Signal Processing， 2020： 8614-8618.

［45］CHEN Y X， CHI Y J. Spectral compressed sensing via structured matrix completion［C］∥Proc.of the 30th International Conference on Machine Learning， 2013： 414-422.

［46］CHEN J L， ZHANG T X， LI J Q， et al. Joint sensor failure detection and corrupted covariance matrix recovery in bistatic MIMO radar with impaired arrays［J］. IEEE Sensors Journal， 2019， 19（14）： 5834-5842.

［47］JIANG H M， CHEN J， DONG G M， et al. Study on Hankel matrix-based SVD and its application in rolling element bearing fault diagnosis［J］. Mechanical Systems and Signal Processing， 2015， 52-53： 338-359.

［48］FAZEL M， PONG T K， SUN D F， et al. Hankel matrix rank minimization with applications to system identification and realization［J］. SIAM Journal on Matrix Analysis and Applications， 2013， 34（3）： 946-977.

作者簡介

張明龍（2000—），男，碩士研究生，主要研究方向為陣列信號處理。

吳雨林（2000—），男，碩士研究生，主要研究方向為毫米波雷達信號處理。

魏文強（1997—），男，博士研究生，主要研究方向為陣列信號處理、最優化理論和算法。

沈園杰（1999—）， 男，碩士研究生，主要研究方向為MIMO毫米波雷達、信號處理。

郭世盛（1991—），男，副研究員，博士，主要研究方向為城市環境多徑探測、雷達信號處理。

崔國龍（1982—），男，教授，博士，主要研究方向為城市環境目標探測、陣列信號處理。