CN110579737B - Sparse array-based MIMO radar broadband DOA calculation method in clutter environment - Google Patents

Abstract

The invention discloses a sparse array-based MIMO radar broadband DOA calculation method in a clutter environment, and belongs to the field of signal processing. Specifically WCSAB, in the method, clutter interference is suppressed by Capon beam forming, and then target DOA is estimated by a CS method by jointly utilizing different narrowband signal information. Considering that the DOA estimation performance is not only related to the beam forming weight value but also related to the sparse array structure, the invention provides a joint optimization problem of the beam forming weight value and the sparse array, and provides a simple algorithm for solving the optimization problem. The method provided by the invention can improve the performance of target DOA estimation in a clutter environment, including high resolution and low sidelobe, the sparse array reduces the system cost and complexity, the Bayesian Mean Square Error (BMSE) of the target DOA estimation is taken as a performance evaluation index, and the sparse array structure designed by the algorithm is similar to the optimal sparse array performance obtained by an exhaustive method and has better performance than a nested array and a co-prime array.

Description

Translated fromChinese

一种杂波环境中基于稀疏阵列的MIMO雷达宽带DOA计算方法A sparse array-based MIMO radar wideband DOA calculation method in clutter environment

技术领域technical field

本发明属于信号处理领域，它特别涉及杂波环境中基于稀疏阵列的MIMO雷达宽带DOA估计问题。The invention belongs to the field of signal processing, and particularly relates to the wideband DOA estimation problem of MIMO radar based on sparse array in clutter environment.

背景技术Background technique

MIMO(Multiple Input Multiple Output)雷达是利用多个发射天线同步地发射信号，同时使用多个接收天线接收回波信号，并集中处理的一种新型雷达系统。相比于传统的相控阵雷达，MIMO雷达具有明显优势，如更高的分辨率，更好的目标检测、定位和跟踪性能，更好的目标参数估计和识别能力。DOA估计研究是阵列信号处理中的一项重要内容，其应用涉及雷达、通信、声呐、射电天文、勘测、地震以及生物医学等领域。目前已有多种经典的DOA估计方法，如多信号分类(Multiple Signal Classification,MUSIC)，基于旋转不变技术的信号参数估计(Estimation of Signal Parameters via Rotational invariancetechniques,ESPRIT)等。近年来，压缩感知(Compressive Sensing,CS)理论得到了国内外学者的广泛关注，相对于传统的方法，基于CS的MIMO雷达DOA估计在较少采样数据和低信噪比条件下拥有更好的估计性能。MIMO (Multiple Input Multiple Output) radar is a new type of radar system that uses multiple transmitting antennas to transmit signals synchronously, and simultaneously uses multiple receiving antennas to receive echo signals and process them centrally. Compared with traditional phased array radar, MIMO radar has obvious advantages, such as higher resolution, better target detection, localization and tracking performance, and better target parameter estimation and identification ability. DOA estimation research is an important content in array signal processing, and its applications involve radar, communication, sonar, radio astronomy, surveying, seismic, and biomedical fields. At present, there are many classical DOA estimation methods, such as Multiple Signal Classification (MUSIC), Estimation of Signal Parameters via Rotational invariance techniques (ESPRIT) and so on. In recent years, Compressive Sensing (CS) theory has received extensive attention from scholars at home and abroad. Compared with traditional methods, CS-based MIMO radar DOA estimation has better performance under the condition of less sampled data and low signal-to-noise ratio. Estimated performance.

根据传统的阵列信号处理理论，为了保证DOA估计的唯一性，阵列中相邻阵元间距要小于或等于入射信号半波长，满足这一条件的阵列称为满阵。阵列的空间分辨率与阵列孔径有关，提高分辨率需要增大阵列孔径，在满阵中，这意味着需要更多的天线。然而，由于实际中软硬件资源的约束，天线的数目通常是受限的。为了增大阵列孔径，同时不增加天线个数，稀疏阵列引起了广泛的关注。当目标在观测空间是稀疏时，稀疏阵列能够准确的估计出目标DOA。然而，在杂波环境中，目标在观测空间的稀疏性会遭到破坏，从而导致DOA估计性能下降。According to the traditional array signal processing theory, in order to ensure the uniqueness of DOA estimation, the distance between adjacent array elements in the array should be less than or equal to the half wavelength of the incident signal. An array that satisfies this condition is called a full array. The spatial resolution of the array is related to the array aperture. Increasing the resolution requires increasing the array aperture. In a full array, this means that more antennas are required. However, due to the constraints of hardware and software resources in practice, the number of antennas is usually limited. In order to increase the array aperture without increasing the number of antennas, sparse arrays have attracted extensive attention. When the target is sparse in the observation space, the sparse array can accurately estimate the target DOA. However, in a cluttered environment, the sparsity of the target in the observation space is destroyed, resulting in the degradation of DOA estimation performance.

考虑到宽带信号具有信息量大、抗干扰能力强、分辨率高等优点，比较有代表性的宽带DOA估计方法有非相干信号子空间(ISSM)和相干信号子空间(CSSM)等。ISSM将宽带信号分成频带上的多个窄带信号，然后分别对每个窄带信号进行处理，最后对所有窄带信号处理的结果求平均得到最终的估计结果，这种方法计算量大，且在低信噪比时估计性能差。CSSM是通过聚焦将不同频率窄带信号的协方差矩阵变换到一个参考频率，然后再用窄带估计方法得到最终结果。但是CSSM需要对目标DOA预估，且性能受预估精度的影响很大。Considering that wideband signals have the advantages of large amount of information, strong anti-interference ability and high resolution, the more representative broadband DOA estimation methods are incoherent signal subspace (ISSM) and coherent signal subspace (CSSM). ISSM divides the wideband signal into multiple narrowband signals on the frequency band, then processes each narrowband signal separately, and finally averages the results of all narrowband signal processing to obtain the final estimation result. The estimation performance is poor at the noise ratio. CSSM transforms the covariance matrix of narrowband signals of different frequencies to a reference frequency by focusing, and then uses the narrowband estimation method to obtain the final result. However, CSSM needs to estimate the target DOA, and the performance is greatly affected by the prediction accuracy.

发明内容SUMMARY OF THE INVENTION

本发明提供了一种杂波环境中基于稀疏阵列的MIMO雷达宽带DOA估计方法，具体为WCSAB(wideband compressive sensing after beamforming)，在该方法中，用Capon波束形成来抑制杂波干扰，然后联合利用不同窄带信号信息，用CS方法估计目标DOA。考虑到DOA估计性能不仅与波束形成权重值有关，还与稀疏阵列结构有关，本发明提出了波束形成权重值和稀疏阵列的联合优化问题，并给出了一种简单的算法求解这一优化问题。The present invention provides a sparse array-based MIMO radar wideband DOA estimation method in a clutter environment, specifically WCSAB (wideband compressive sensing after beamforming). In this method, Capon beamforming is used to suppress clutter interference, and then combined use Different narrowband signal information, the CS method is used to estimate the target DOA. Considering that the DOA estimation performance is not only related to the beamforming weight value, but also to the sparse array structure, the present invention proposes a joint optimization problem of the beamforming weight value and the sparse array, and provides a simple algorithm to solve this optimization problem .

本发明技术方案为一种杂波环境中基于稀疏阵列的MIMO雷达宽带DOA计算方法，该方法包括：The technical solution of the present invention is a sparse array-based MIMO radar wideband DOA calculation method in a clutter environment, the method comprising:

步骤1：设发射天线位置确定，放置接收天线的可行域为[0,D_r]，为了简化分析，将可行域以间隔Δ_r离散化为N_r个格点，且有N个接收天线放置在其中一些格点上，N＜＜N_r；Step 1: Assume that the position of the transmitting antenna is determined, and the feasible domain for placing the receiving antenna is [0, D_r ]. In order to simplify the analysis, the feasible domain is discretized into N_r lattice points with an interval Δ_r , and there are N receiving antennas placed. On some of the lattice points, N<<N_r ;

步骤2：建立MIMO雷达回波信号模型，得到回波信号时域采样数据

n＝1,...,N_r和p＝1,...,L，其中p表示时域快拍，L为快拍数；Step 2: Establish a MIMO radar echo signal model to obtain time-domain sampling data of the echo signal

n=1,...,N_r and p=1,...,L, where p represents the snapshot in the time domain, and L is the number of snapshots;

步骤3：对接收信号

进行L点离散傅里叶变换得到频域数据，即

Step 3: On the Received Signal

Perform L-point discrete Fourier transform to obtain frequency domain data, namely

并将N_r个格点的数据表示成矢量形式，即y[l]＝[y₁[l],...,y_Nr[l]]^T，其中p＝1,...,L，l＝1,...,L；And represent the data of N_r lattice points in vector form, that is, y[l]=[y₁ [l],...,y_Nr [l]]^T , where p=1,...,L, l=1,...,L;

步骤4：将目标角度观测区域离散化为G个格点θ₁,...,θ_G，K＜＜G，其中K表示目标个数，将信号模型表示成稀疏形式：Step 4: Discretize the target angle observation area into G lattice points θ₁ ,...,θ_G , K<<G, where K represents the number of targets, and the signal model is expressed in a sparse form:

y[l]＝Φ[l]x+c[l]+u[l]y[l]=Φ[l]x+c[l]+u[l]

其中

这里a_r(θ,f_l)表示接收导向矢量，a_t(θ,f_l)表示发射导向矢量，s[l]表示频域发射信号，x＝[x₁,...,x_G]^T是K稀疏的，也就是x只有K个非零元素，且非零元素的值和位置为目标反射系数和DOA，c[l]表示杂波，u[l]表示噪声；in

Here ar (θ,f_l ) represents the receive steering_vector , at (θ,_fl ) represents the transmit steering vector, s[l] represents the frequency domain transmit signal, x=[x₁ ,...,x_G ]^T is K sparse, that is, x has only K non-zero elements, and the value and position of the non-zero elements are the target reflection coefficient and DOA, c[l] represents clutter, and u[l] represents noise;

步骤5：将波束形成权重矢量w_g,l作用到y[l]上得到波束形成输出结果：Step 5: Apply the beamforming weight vector w_g,l to y[l] to get the beamforming output result:

将r_g,l(g＝1,...,G，l＝1,...,L)表示成矢量：

Represent r_g,l (g=1,...,G,l=1,...,L) as a vector:

r＝[r_1,1,...,r_G,1,...,r_1,L,...,r_G,L]^Tr=[r_1,1 ,...,r_G,1 ,...,r_1,L ,...,r_G,L ]^T

＝W_rΦx+W_rc+W_ru=W_r Φx+W_r c+W_r u

其中权值矩阵W_r＝Diag{W₁,...,W_l,...,W_L}是一个块对角矩阵；Wherein the weight matrix W_r =Diag{W₁ ,...,W_l ,...,W_L } is a block diagonal matrix;

且有

Φ＝[Φ^T[1],...,Φ^T[L]]^T，c＝[c^T[1],...,c^T[L]]^T表示杂波，u＝[u^T[1],...,u^T[L]]^T表示噪声；and have

Φ=[Φ^T [1],...,Φ^T [L]]^T , c=[c^T [1],...,c^T [L]]^T represents clutter, u=[u^T [1],...,u^T [L]]^T represents noise;

步骤6：基于CS理论，通过基寻踪去燥来重构稀疏向量x；Step 6: Based on CS theory, the sparse vector x is reconstructed by base pursuit to remove dryness;

其中η≥0是正则化参数；where η≥0 is the regularization parameter;

步骤7：对步骤6得到的解

中的元素值按照从大到小排序，排序后各个元素相应的格点表示为{θ₍₁₎,...,θ_(G)}，那么DOA估计结果可以表示为；Step 7: Solution toStep 6

The element values in are sorted from large to small. After sorting, the corresponding grid points of each element are expressed as {θ₍₁₎ ,...,θ_(G) }, then the DOA estimation result can be expressed as;

步骤8：基于最小化贝叶斯均方误差

求解最优W_r，建立以下优化问题Step 8: Minimize Bayesian mean squared error based on

To solve the optimal W_r , the following optimization problem is established

s.t.W_r＝Diag{W₁,...,W_L}stW_r =Diag{W₁ ,...,W_L }

||W||₀＝N||W||₀ = N

W＝[w_1,1,...,w_1,L,...,w_G,1,...,w_G,L]W=[w_1,1 ,...,w_1,L ,...,w_G,1 ,...,w_G,L ]

其中，真实目标的DOA矢量θ_T是随机的，

表示对θ_T求期望，

表示θ_T确定时，DOA估计的均方误差，w_g,l表示权重矢量，其中g＝1,...,G，l＝1,...,L；Among them, the DOA vector θ_T of the real target is random,

represents the expectation for θ_T ,

Represents the mean square error of DOA estimation when θ_T is determined, w_g,l represents the weight vector, where g=1,...,G, l=1,...,L;

步骤9：优化求解步骤8提出的问题，得到最优的W_r。Step 9: Optimally solve the problem raised inStep 8 to obtain the optimal W_r .

进一步的，所述步骤9的具体方法为：Further, the specific method of thestep 9 is:

步骤1：初始化：迭代次数j＝1，

根据公式

计算波束形成权重值

其中

R_c(f_l)为杂波c[l]的协方差矩阵；在每次迭代过程中，随机产生一组格点选择矢量{z₁,...,z_α}，对于给定的z；Step 1: Initialization: The number of iterations j=1,

According to the formula

Calculate the beamforming weight value

in

R_c (f_l ) is the covariance matrix of the clutter c[l]; in each iteration, a set of lattice point selection vectors {z₁ ,...,z_α } are randomly generated, for a given z ;

步骤2：重复步骤3到步骤6的迭代过程：Step 2: Repeat the iterative process fromsteps 3 to 6:

步骤3：随机产生一组格点选择矢量{z₁,...,z_α}；Step 3: Randomly generate a set of lattice point selection vectors {z₁ ,...,z_α };

步骤4：根据公式w_g,l＝z⊙ξ_g,l计算

并构成

Step 4: Calculate according to the formula w_g,l =z⊙ξg_,l

and constitute

根据公式

计算r_g,l，将r_g,l表示成矢量r；According to the formula

Calculate r_g,l and represent r_g,l as a vector r;

将

和r代入公式

得到x重构结果

和目标DOA估计结果

Will

and r into the formula

get x reconstruction result

and target DOA estimation results

根据公式

得到BMSE

According to the formula

get BMSE

步骤5：基于最小BMSE得到

Step 5: get based on the minimum BMSE

步骤6：根据

基于

更新

中相应的权重值得到

并令j＝j+1；其中

为Dc[l]的协方差矩阵，

Step 6: According to

based on

renew

The corresponding weight values in the

and let j=j+1; where

is the covariance matrix of Dc[l],

步骤7：当

时，迭代停止，输出最优天线选择

e₀为事先设定的阈值。Step 7: When

When , the iteration stops and the optimal antenna selection is output

e₀ is a preset threshold.

本发明提出的方法可以提高杂波环境中目标DOA估计的性能，包括高分辨率和低旁瓣，稀疏阵列降低了系统成本和复杂性，以目标DOA估计的贝叶斯均方误差(BMSE)为性能评价指标，通过上述算法设计的稀疏阵列结构与穷举法得到的最优稀疏阵列性能相近，且比嵌套阵和互质阵性能更优。The method proposed in the present invention can improve the performance of target DOA estimation in clutter environment, including high resolution and low side lobes, sparse array reduces system cost and complexity, Bayesian mean square error (BMSE) for target DOA estimation As a performance evaluation index, the sparse array structure designed by the above algorithm has similar performance to the optimal sparse array obtained by the exhaustive method, and is better than the nested array and the coprime array.

附图说明Description of drawings

图1给出了所有可能稀疏阵列结构下的BMSE以升序排列的结果，为了对比，图1也给出了嵌套阵(nested array)和互质阵(co-prime array)的结果。Figure 1 presents the results of BMSE in ascending order for all possible sparse array structures. For comparison, Figure 1 also presents the results for nested arrays and co-prime arrays.

图2(a)给出了基于最小BMSE条件下最优的稀疏阵列结构，图2(b)给出了根据本发明所提出的算法得到的稀疏阵列结构。Fig. 2(a) shows the optimal sparse array structure under the condition of minimum BMSE, and Fig. 2(b) shows the sparse array structure obtained by the algorithm proposed in the present invention.

图3给出了使用WCSAB方法时，不同稀疏阵列结构的DOA估计结果。Figure 3 presents the DOA estimation results for different sparse array structures when using the WCSAB method.

图4考虑单目标情况，分别用WCSAB和WCT(wideband Capon technique)方法时，不同阵列结构的DOA估计结果。Figure 4 considers the single target case, when using the WCSAB and WCT (wideband Capon technique) methods respectively, the DOA estimation results of different array structures.

图5为考虑双目标情况，分别用WCSAB和WCT(wideband Capon technique)方法时，不同阵列结构的DOA估计结果。Figure 5 shows the DOA estimation results of different array structures when WCSAB and WCT (wideband Capon technique) methods are used, respectively, considering dual targets.

具体实施方式Detailed ways

为了方便描述，首先进行如下定义：For the convenience of description, the following definitions are first made:

黑体大写字母表示矩阵，黑体小写字母表示矢量，(·)^*为共轭，(·)^T为转置，(·)^H为共轭转置，||x||₀和||x||₁分别表示向量x的l₀范数和l₁范数，||W||₀表示矩阵W非零行的个数，Diag{·}表示块对角矩阵，diag_r{·}表示去掉零行之后的对角矩阵，

表示相对于θ的期望，I_N为N阶的单位阵，1为全1矢量，符号⊙表示哈达玛积。Bold uppercase letters indicate matrices, bold lowercase letters indicate vectors, ( )^* is conjugate, ( )^T is transpose, ( )^H is conjugate transpose, ||x||₀ and ||x||₁ represents the l₀ norm and l₁ norm of the vector x, respectively, ||W||₀ represents the number of non-zero rows in the matrix W, Diag{·} represents a block diagonal matrix, and diag_r {·} means remove zeros the diagonal matrix after the row,

Represents the expectation relative to θ, I_N is a unit matrix of order N, 1 is an all-one vector, and the symbol ⊙ represents the Hadamard product.

考虑一个共置MIMO雷达系统，发射天线和接收天线都放置在二维笛卡尔坐标系的横轴上。假设有M个发射天线，且在横轴上的位置已知，为d_t,m(m＝1,...,M)。假设放置接收天线的可行域为[0,D_r]，为了简化分析，将可行域以间隔Δ_r离散化为N_r个格点，接收天线放置在这些格点上。由于天线个数的约束，假设雷达系统只有N(N＜＜N_r)个可用的接收天线。令

表示第m个发射天线发射的宽带信号，频率范围为[-B_m/2,B_m/2]，其中p表示时域快拍，T_s表示采样周期，L表示快拍数。假设K个远场点目标的DOA为θ_T,k(k＝1,...,K)，那么在第n个格点接收到的信号为Consider a co-located MIMO radar system with both transmit and receive antennas placed on the horizontal axis of a two-dimensional Cartesian coordinate system. It is assumed that there are M transmit antennas, and the positions on the horizontal axis are known as d_t,m (m=1,...,M). Assuming that the feasible region where the receiving antenna is placed is [0, D_r ], in order to simplify the analysis, the feasible region is discretized into N_r lattice points at the interval Δ_r , and the receiving antenna is placed on these lattice points. Due to the constraint of the number of antennas, it is assumed that the radar system has only N (N<<N_r ) available receive antennas. make

Indicates the broadband signal transmitted by the mth transmitting antenna, the frequency range is [-B_m /2, B_m /2], where p represents the time domain snapshot, T_s represents the sampling period, and L represents the number of snapshots. Assuming that the DOA of the K far-field point targets is θ_T,k (k=1,...,K), then the signal received at the nth grid point is

其中f_c表示载频，β_k表示第k个目标的反射系数，且假设是确定未知的。令第一个发射天线和第一个格点作为参考，那么τ_Tt,m,k＝(d_t,m-d_t,1)sinθ_T,k/c表示信号从第m个发射天线到第k个目标时，相对于参考阵元的时延，τ_Tr,n,k＝(n-1)Δ_rsinθ_T,k表示信号从第k个目标到第n个格点时，相对于第一个格点的时延。Q表示杂波散射体的个数，γ_q(q＝1,...,Q)表示杂波散射体的反射系数，并且假设它们之间是独立同分布(iid)的高斯随机变量。τ_Ct,m,q＝(d_t,m-d_t,1)sinθ_C,q/c表示信号从第m个发射天线到第q个杂波散射体时，相对于参考阵元的时延，τ_Cr,n,q＝(n-1)Δ_rsinθ_C,q表示信号从第q个杂波散射体到第n个格点时，相对于第一个格点的时延，θ_C,q表示第q个杂波散射体相对于阵列的方向。

是方差为σ²的高斯白噪声。where f_c is the carrier frequency, β_k is the reflection coefficient of the k-th target, and the assumption is that it is deterministically unknown. Let the first transmitting antenna and the first lattice point be used as a reference, then τ_Tt,m,k =(d_t,m -d_t,1 )sinθ_T,k /c means that the signal goes from the mth transmitting antenna to the When there are k targets, relative to the time delay of the reference array element, τ_Tr,n,k =(n-1)Δ_r sinθ_T,k indicates that when the signal goes from the kth target to the nth grid point, relative to the One grid point delay. Q represents the number of clutter scatterers, γ_q (q=1, . τ_Ct,m,q =(d_t,m -d_t,1 )sinθ_C,q /c represents the time delay relative to the reference array element when the signal travels from the mth transmitting antenna to the qth clutter scatterer , τ_Cr,n,q =(n-1)Δ_r sinθ_C,q represents the time delay relative to the first grid point when the signal goes from the qth clutter scatterer to the nth grid point, θ_{C ,q} denotes the direction of the qth clutter scatterer relative to the array.

is white Gaussian noise with variance^σ2 .

通过对时域离散信号进行L点离散傅里叶变换(DFT)，可得到在频率点f_l＝lf_s(l＝1,...,L)的频域数据，其中f_s为频率采样间隔，f_l∈[-B/2,B/2]且

信号在频率f_l处的DFT结果为By performing L-point discrete Fourier transform (DFT) on the time-domain discrete signal, the frequency-domain data at the frequency point f_l =lf_s (l=1,...,L) can be obtained, where f_s is the frequency sampling interval, f_l ∈ [-B/2,B/2] and

The DFT result of the signal at frequency f_l is

其中s_m[l]和u_n[l]分别表示发射信号

和噪声

的DFT。令

和

分别表示在角度θ、频率f_l处的接收导向矢量和发射导向矢量。将N_r个格点接收到的信号表示成矢量where s_m [l] and u_n [l] represent the transmitted signal, respectively

and noise

DFT. make

and

denote the receive steering vector and transmit steering vector at angle θ and frequency_fl , respectively. Represent the signals received at N_r lattice points as a vector

其中

in

在CS框架下，为了估计K个目标的DOA θ_T,k(k＝1,...,K)，将目标角度观测区域离散化为G(K＜＜G)个格点θ₁,...,θ_G，假设离散误差可以忽略，即目标正好落在格点上。那么(3)式可以表示为Under the CS framework, in order to estimate DOA θ_T,k (k=1,...,K) of K targets, the target angle observation area is discretized into G (K<<G) lattice points θ₁ ,. ..,θ_G , assuming that the discrete error is negligible, that is, the target falls exactly on the grid point. Then (3) can be expressed as

y[l]＝Φ[l]x+c[l]+u[l] (4)y[l]=Φ[l]x+c[l]+u[l] (4)

其中

矢量x＝[x₁,...,x_G]^T是K稀疏的，也就是x只有K个非零元素，且非零元素的值和位置为目标反射系数和DOA，可以表示为in

The vector x=[x₁ ,...,x_G ]^T is K sparse, that is, x has only K non-zero elements, and the value and position of the non-zero elements are the target reflection coefficient and DOA, which can be expressed as

CS理论利用x的稀疏性来估计目标DOA，然而，这种稀疏性在杂波环境中会遭到破坏，从而降低DOA估计的性能。为了抑制杂波的干扰，在接收端采用波束形成的方法。令

表示在方向θ_g、频率f_l处的波束形成权重矢量，且非零元素的位置表示选择放置天线的格点。由于只有N个可用的接收天线，因此要求权重矢量满足||w_g,l||₀＝N。波束形成的输出由下式给出CS theory exploits the sparsity of x to estimate the target DOA, however, this sparsity is corrupted in a cluttered environment, degrading the performance of DOA estimation. In order to suppress the interference of clutter, the beamforming method is adopted at the receiving end. make

represents the beamforming weight vector at direction θ_g , frequency_fl , and the position of the non-zero element represents the grid point where the antenna is chosen to be placed. Since there are only N available receive antennas, the weight vector is required to satisfy ||w_g,l ||₀ =N. The output of beamforming is given by

将r_g,l(g＝1,...,G和l＝1,...,L)表示成一个GL×1的矢量Represent r_g,l (g=1,...,G and l=1,...,L) as a GL×1 vector

其中Φ＝[Φ^T[1],...,Φ^T[L]]^T，c＝[c^T[1],...,c^T[L]]^T，u＝[u^T[1],...,u^T[L]]^T，W_r＝Diag{W₁,...,W_L}，

根据(7)式，DOA估计问题可以转化为稀疏信号重构问题，基于CS理论，K稀疏的矢量x可以通过基寻踪去燥(BPDN)来重构where Φ=[Φ^T [1],...,Φ^T [L]]^T , c=[c^T [1],...,c^T [L]]^T , u=[u^T [1 ],...,u^T [L]]^T , W_r =Diag{W₁ ,...,W_L },

According to Eq. (7), the DOA estimation problem can be transformed into a sparse signal reconstruction problem. Based on CS theory, the K-sparse vector x can be reconstructed by basis pursuit denoising (BPDN).

其中η≥0是正则化参数，对于(8)式这个优化问题，可以使用CVX工具包求解。令

表示上式的解，

为目标DOA的估计结果。考虑目标DOA矢量θ_T＝[θ_T,1,...,θ_T,K]^T是随机的情况，那么平均估计性能可以由贝叶斯均方误差(BMSE)给出where η≥0 is the regularization parameter. For the optimization problem of (8), the CVX toolkit can be used to solve it. make

represents the solution of the above equation,

is the estimated result of the target DOA. Considering the case where the target DOA vector θ_T = [θ_T,1 ,...,θ_T,K ]^T is random, then the average estimation performance can be given by the Bayesian mean squared error (BMSE)

由(9)式可知，DOA估计的性能与矩阵W_r有关，为了使性能最优，给出以下优化问题It can be seen from equation (9) that the performance of DOA estimation is related to the matrix W_r . In order to optimize the performance, the following optimization problem is given

(10)式中最后两个约束是为了保证对于不同的g和l，w_g,l中非零元素的位置是相同的。由于w_g,l中非零元素的位置表示选中相应的格点放置天线，因此(10)式是一个权重值和稀疏阵列结构的联合优化问题。The last two constraints in (10) are to ensure that for different g and l, the positions of non-zero elements in w_g,l are the same. Since the position of the non-zero element in w_g,l indicates that the corresponding grid point is selected to place the antenna, equation (10) is a joint optimization problem of weight value and sparse array structure.

考虑到(10)式是一个NP-hard问题，提出一种简单的算法来求解该优化问题。该算法首先在权值给定时优化稀疏阵列结构，然后再去更新权值用于下一次迭代。首先解释如何在权值给定时优化稀疏阵列结构。首先，定义一个格点选择矢量

其中z_n∈{0,1}，只有元素为1时表示选择相应的格点放置天线，由于只有N个可用接收天线，要求||z||₀＝N。第一次迭代的权值由满阵情况Capon波束形成给出Considering that Equation (10) is an NP-hard problem, a simple algorithm is proposed to solve the optimization problem. The algorithm first optimizes the sparse array structure when the weights are given, and then updates the weights for the next iteration. First explain how to optimize the sparse array structure when the weights are given. First, define a grid selection vector

Among them, z_n ∈ {0,1}, only when the element is 1 indicates that the corresponding grid point is selected to place the antenna. Since there are only N available receiving antennas, ||z||₀ =N is required. The weights for the first iteration are given by the full array case Capon beamforming

其中

R_c(f_l)为杂波c[l]的协方差矩阵。在每次迭代过程中，随机产生一组格点选择矢量{z₁,...,z_α}，对于给定的z，有in

R_c (f_l ) is the covariance matrix of the clutter c[l]. During each iteration, a set of lattice selection vectors {z₁ ,...,z_α } are randomly generated. For a given z, we have

w_g,l＝z⊙ξ_g,l (12)w_g,l = z⊙ξ_g,l (12)

对于不同的z，可以得到不同的w_g,l和W_r(W_r由w_g,l构成)，由(9)式可知，BMSE与W_r有关，因此可以发现，BMSE与z也相关，表示为e(z)。基于最小BMSE，可以得到最优的格点选择矢量z_op。For different z, different w_{g, l} and W_r can be obtained (W_r is composed of w_{g, l} ). From formula (9), it can be seen that BMSE is related to W_r , so it can be found that BMSE is also related to z, Denoted as e(z). Based on the minimum BMSE, the optimal lattice point selection vector z_op can be obtained.

下面基于z_op更新ξ_g,l中相应的权重值，令

则ξ_g,l中由z_op选中的元素值通过下式更新Next, update the corresponding weight values in ξ_g,l based on z_op , let

Then the element value selected by z_op in ξ_g,l is updated by the following formula

其中

为Dc[l]的协方差矩阵，

当BMSE e(z_op)小于某个阈值e₀时，迭代停止。详细的算法由表1给出。in

is the covariance matrix of Dc[l],

The iteration stops when BMSE e(z_op ) is less than a certain threshold e₀ . The detailed algorithm is given in Table 1.

表1求解优化问题的迭代算法Table 1 Iterative Algorithms for Solving Optimization Problems

为了抑制杂波的干扰，在接收端采用波束形成的方法。

表示在方向θ_g、频率f_l处的波束形成权重矢量，且非零元素的位置表示选择放置天线的格点，由于只有N个可用的接收天线，因此要求权重矢量满足||w_g,l||₀＝N。对于不同的g和l，天线的阵列结构应是相同的，也就是w_g,l中非零元素的位置是相同的。为了表示这一约束，构造矩阵In order to suppress the interference of clutter, the beamforming method is adopted at the receiving end.

Represents the beamforming weight vector at the direction θ_g and frequency f_l , and the position of the non-zero element indicates the grid point where the antenna is placed. Since there are only N available receiving antennas, the weight vector is required to satisfy ||w_g,l ||₀ =N. For different g and l, the array structure of the antenna should be the same, that is, the positions of non-zero elements in w_g,l are the same. To represent this constraint, construct the matrix

W＝[w_1,1,...,w_1,L,...,w_G,1,...,w_G,L] (14)W=[w_1,1 ,...,w_1,L ,...,w_G,1 ,...,w_G,L ] (14)

且满足||W||₀＝N，即矩阵非零行的个数为N，通过这个约束，即可满足对于不同的g和l，w_g,l中非零元素的位置相同。将权重矢量w_g,l作用于接收信号y[l]，根据(6)式得到波束形成输出r_g,l，将其表示成G×1的矢量And satisfy ||W||₀ =N, that is, the number of non-zero rows in the matrix is N, through this constraint, it can be satisfied that for different g and l, the positions of non-zero elements in w_{g, l} are the same. The weight vector w_g,l is applied to the received signal y[l], and the beamforming output r_g,l is obtained according to equation (6), which is expressed as a vector of G×1

由上式可以发现，对于不同的频率f_l，向量x稀疏性一样。为了联合利用不同频率的信号信息，将r_l(l＝1,...,L)表示成GL×1的矢量，即(7)式。通过将DOA估计问题可以转化为稀疏信号重构问题，根据(8)式，可以得到稀疏向量x的重构结果

则

中最大的K个元素的位置为目标DOA的估计结果，表示为

对

中的元素值按照从大到小排序，排序后各个元素相应的格点表示为{θ₍₁₎,...,θ_(G)}，那么DOA估计结果可以表示为It can be found from the above formula that for different frequencies f_l , the vector x has the same sparsity. In order to jointly utilize the signal information of different frequencies, r_l (l=1, . . . , L) is expressed as a GL×1 vector, that is, formula (7). By transforming the DOA estimation problem into a sparse signal reconstruction problem, according to equation (8), the reconstruction result of the sparse vector x can be obtained

but

The position of the largest K elements in is the estimation result of the target DOA, which is expressed as

right

The element values in are sorted from large to small. After sorting, the corresponding grid points of each element are expressed as {θ₍₁₎ ,...,θ_(G) }, then the DOA estimation result can be expressed as

满阵波束形成权重矢量ξ_g,l可以使方向θ_g、频率f_l处的信号无失真通过，同时抑制其他方向的干扰和噪声，表示为The full array beamforming weight vector ξ_g,l can make the signal at the direction θ_g and frequency f_l pass without distortion, while suppressing the interference and noise in other directions, which is expressed as

上式的最优解即为(11)式。The optimal solution of the above equation is equation (11).

关于杂波环境中基于稀疏阵列的MIMO雷达宽带DOA估计，给出了两个仿真实例，参数设置如下：假设D_r＝11λ/2，其中λ表示信号最高频率对应的波长。将可放置接收天线的可行域以Δ_r＝λ/2为间隔离散为12个格点。假设MIMO雷达系统可用的发射和接收天线个数为M＝N＝6，且发射端阵列结构确定已知。Regarding the wideband_DOA estimation of MIMO radar based on sparse array in clutter environment, two simulation examples are given. The feasible region where the receiving antenna can be placed is discretized into 12 lattice points at the interval of Δ_r =λ/2. It is assumed that the number of available transmit and receive antennas in the MIMO radar system is M=N=6, and the array structure of the transmitting end is determined and known.

为了简化分析，假设发射信号带宽相同，即B_m＝200MHz(m＝1,...,M)，载频为1GHz。In order to simplify the analysis, it is assumed that the transmitted signal bandwidth is the same, that is, B_m =200MHz (m=1, . . . , M), and the carrier frequency is 1GHz.

将目标角度观测域离散为41个格点-20°,-19°,...,20°。The target angle observation domain is discretized into 41 grid points -20°,-19°,...,20°.

杂波由250个散射体构成，分布的角度为-90°,-90°+180°/250,...,90°。The clutter consists of 250 scatterers distributed at angles of -90°, -90°+180°/250,...,90°.

定义信噪比

和信杂比

不失一般性，假设目标反射系数为1，SNR和SCR分别设置为-5dB和-30dB。Define the signal-to-noise ratio

and signal-to-noise ratio

Without loss of generality, it is assumed that the target reflection coefficient is 1, and the SNR and SCR are set to -5dB and -30dB, respectively.

在仿真1中，假设目标在角度观测域离散化后的格点上是均匀随机分布的。为了保证阵列孔径不改变，令

那么共有

种不同的稀疏阵列结构。图1给出了所有可能稀疏阵列结构下的BMSE以升序排列的结果，菱形表示最小BMSE条件下最优的稀疏阵列结构，正方形表示根据表1所给出的算法得到的稀疏阵列结构，可以看到这两种结构的性能相近，这两种稀疏阵列的具体结构由图2给出。为了对比，图1也给出了嵌套阵(nested array)和互质阵(co-prime array)的结果，分别由星号和圆形标出，可以看到有多种稀疏阵列结构的性能都比嵌套阵和互质阵更好。假设只有一个目标，DOA为-14°，图3给出了使用WCSAB方法时，上述四种稀疏阵列结构的DOA估计结果。从图中可以看到，最优稀疏阵列和由算法得到的稀疏阵列可以准确估计目标DOA，而嵌套阵和互质阵估计有误差。Insimulation 1, it is assumed that the targets are uniformly and randomly distributed on the grid points after the discretization of the angle observation domain. In order to ensure that the array aperture does not change, let

then share

different sparse array structures. Figure 1 shows the results of BMSEs under all possible sparse array structures arranged in ascending order. The diamonds represent the optimal sparse array structures under the condition of minimum BMSE, and the squares represent the sparse array structures obtained according to the algorithm given in Table 1. It can be seen that The performance of these two structures is similar, and the specific structures of these two sparse arrays are given in Figure 2. For comparison, Figure 1 also shows the results of nested array and co-prime array, which are marked by asterisks and circles, respectively. It can be seen that there are various performances of sparse array structures. Both are better than nested and coprime matrices. Assuming that there is only one target and the DOA is -14°, Figure 3 shows the DOA estimation results of the above four sparse array structures when using the WCSAB method. It can be seen from the figure that the optimal sparse array and the sparse array obtained by the algorithm can accurately estimate the target DOA, while the estimation of the nested array and the coprime array has errors.

在仿真2中，比较了WCSAB与WCT(wideband Capon technique)这两种宽带DOA估计方法的性能。WCT方法属于ISSM的一种，它对每一个窄带信号用Capon方法得到相应的DOA估计结果，然后对所有的结果求平均得到最终的估计结果。图4考虑单目标情况，目标DOA为10°，图5考虑双目标情况，目标DOA为6°和10°。图中实线表示以λ/2为间隔的满阵结构，虚线表示最优稀疏阵列结构，点划线表示根据算法得到的稀疏阵列结构。从图4和图5可以看出，相比于满阵结构，两种稀疏阵列所对应的主瓣宽度更窄，也就是说稀疏阵列的分辨率更高。但是稀疏阵列会导致较高的旁瓣，并且这一问题在双目标情况下更严重。通过比较，可以发现，WCSAB方法的旁瓣比WCT低。由图5也可以看到，WCSAB方法在两种稀疏阵列中都可以准确的估计出目标DOA，而WCT估计有误差。通过对比，WCSAB方法的性能更好。InSimulation 2, the performance of two wideband DOA estimation methods, WCSAB and WCT (wideband Capon technique), is compared. The WCT method belongs to a kind of ISSM. It uses the Capon method to obtain the corresponding DOA estimation result for each narrowband signal, and then averages all the results to obtain the final estimation result. Figure 4 considers a single target case with a target DOA of 10°, and Figure 5 considers a dual target case with a target DOA of 6° and 10°. The solid line in the figure represents the full array structure with λ/2 as the interval, the dotted line represents the optimal sparse array structure, and the dot-dash line represents the sparse array structure obtained according to the algorithm. It can be seen from Fig. 4 and Fig. 5 that, compared with the full array structure, the width of the main lobe corresponding to the two sparse arrays is narrower, that is to say, the resolution of the sparse array is higher. But sparse arrays lead to higher sidelobes, and this problem is more severe in the dual-objective case. By comparison, it can be found that the side lobe of the WCSAB method is lower than that of the WCT. It can also be seen from Figure 5 that the WCSAB method can accurately estimate the target DOA in both sparse arrays, while the WCT estimation has errors. By comparison, the performance of the WCSAB method is better.

Claims (2)

Translated fromChinese

1.一种杂波环境中基于稀疏阵列的MIMO雷达宽带DOA计算方法，该方法包括：1. A sparse array-based MIMO radar broadband DOA calculation method in a clutter environment, the method comprising:步骤1：设发射天线位置确定，放置接收天线的可行域为[0,D_r]，为了简化分析，将可行域以间隔Δ_r离散化为N_r个格点，且有N个接收天线放置在其中一些格点上，N＜＜N_r；Step 1: Assume that the position of the transmitting antenna is determined, and the feasible domain for placing the receiving antenna is [0, D_r ]. In order to simplify the analysis, the feasible domain is discretized into N_r lattice points with an interval Δ_r , and there are N receiving antennas placed. On some of the lattice points, N<<N_r ;步骤2：建立MIMO雷达回波信号模型，得到回波信号时域采样数据

n＝1,...,N_r和p＝1,...,L，其中p表示时域快拍，L为快拍数；Step 2: Establish a MIMO radar echo signal model to obtain time-domain sampling data of the echo signal

n=1,...,N_r and p=1,...,L, where p represents the snapshot in the time domain, and L is the number of snapshots;步骤3：对接收信号

进行L点离散傅里叶变换得到频域数据，即

并将N_r个格点的数据表示成矢量形式，即

其中p＝1,...,L，l＝1,...,L；Step 3: On the Received Signal

Perform L-point discrete Fourier transform to obtain frequency domain data, namely

And represent the data of N_r grid points in vector form, that is

where p=1,...,L, l=1,...,L;步骤4：将目标角度观测区域离散化为G个格点θ₁,...,θ_G，K＜＜G，其中K表示目标个数，将信号模型表示成稀疏形式：Step 4: Discretize the target angle observation area into G lattice points θ₁ ,...,θ_G , K<<G, where K represents the number of targets, and the signal model is expressed in a sparse form:y[l]＝Φ[l]x+c[l]+u[l]y[l]=Φ[l]x+c[l]+u[l]其中

这里a_r(θ,f_l)表示接收导向矢量，f_l表示频率点，a_t(θ,f_l)表示发射导向矢量，s[l]表示频域发射信号，x＝[x₁,...,x_G]^T是K稀疏的，也就是x只有K个非零元素，且非零元素的值和位置为目标反射系数和DOA，c[l]表示杂波，u[l]表示噪声；in

Here a_r (θ, f_l ) represents the receive steering vector, f_l represents the frequency point, at (θ,_fl ) represents the transmit steering vector, s[l] represents the frequency domain transmit signal, x=[x₁ ,. ..,x_G ]^T is K sparse, that is, x has only K non-zero elements, and the value and position of the non-zero elements are the target reflection coefficient and DOA, c[l] represents clutter, and u[l] represents noise;步骤5：将波束形成权重矢量w_g,l作用到y[l]上得到波束形成输出结果：Step 5: Apply the beamforming weight vector w_g,l to y[l] to get the beamforming output result:

将r_g,l表示成矢量：

Represent r_g,l as a vector:r＝[r_1,1,...,r_G,1,...,r_1,L,...,r_G,L]^Tr=[r_1,1 ,...,r_G,1 ,...,r_1,L ,...,r_G,L ]^T＝W_rΦx+W_rc+W_ru=W_r Φx+W_r c+W_r u其中权值矩阵W_r＝Diag{W₁,...,W_l,...,W_L}是一个块对角矩阵，g＝1,...,G，l＝1,...,L；where the weight matrix W_r =Diag{W₁ ,...,W_l ,...,W_L } is a block diagonal matrix, g=1,...,G, l=1,... ,L;且有

Φ＝[Φ^T[1],...,Φ^T[L]]^T，c＝[c^T[1],...,c^T[L]]^T表示杂波，u＝[u^T[1],...,u^T[L]]^T表示噪声；and have

Φ=[Φ^T [1],...,Φ^T [L]]^T , c=[c^T [1],...,c^T [L]]^T represents clutter, u=[u^T [1],...,u^T [L]]^T represents noise;步骤6：基于CS理论，通过基寻踪去噪来重构稀疏向量x；Step 6: Based on CS theory, reconstruct the sparse vector x through base pursuit denoising;

其中η≥0是正则化参数；where η≥0 is the regularization parameter;步骤7：对步骤6得到的解

中的元素值按照从大到小排序，排序后各个元素相应的格点表示为{θ₍₁₎,...,θ_(G)}，那么DOA估计结果可以表示为；Step 7: Solution to Step 6

The element values in are sorted from large to small. After sorting, the corresponding grid points of each element are expressed as {θ₍₁₎ ,...,θ_(G) }, then the DOA estimation result can be expressed as;

步骤8：基于最小化贝叶斯均方误差

求解最优W_r，建立以下优化问题Step 8: Minimize Bayesian mean squared error based on

To solve the optimal W_r , the following optimization problem is established

s.t.W_r＝Diag{W₁,...,W_L}stW_r =Diag{W₁ ,...,W_L }

||W||₀＝N||W||₀ = NW＝[w_1,1,...,w_1,L,...,w_G,1,...,w_G,L]W=[w_1,1 ,...,w_1,L ,...,w_G,1 ,...,w_G,L ]其中，真实目标的DOA矢量θ_T是随机的，

表示对θ_T求期望，

表示θ_T确定时，DOA估计的均方误差，w_g,l表示权重矢量，其中g＝1,...,G，l＝1,...,L；Among them, the DOA vector θ_T of the real target is random,

represents the expectation for θ_T ,

Represents the mean square error of DOA estimation when θ_T is determined, w_g,l represents the weight vector, where g=1,...,G, l=1,...,L;步骤9：优化求解步骤8提出的问题，得到最优的W_r。Step 9: Optimally solve the problem raised in Step 8 to obtain the optimal W_r .2.如权利要求1所述的一种杂波环境中基于稀疏阵列的MIMO雷达宽带DOA计算方法，其特征在于所述步骤9的具体方法为：2. in a kind of clutter environment as claimed in claim 1, the MIMO radar wideband DOA calculation method based on sparse array is characterized in that the concrete method of described step 9 is:步骤9.1：初始化：迭代次数j＝1，

根据公式

计算波束形成权重值

其中

R_c(f_l)为杂波c[l]的协方差矩阵；在每次迭代过程中，随机产生一组格点选择矢量{z₁,...,z_α}，对于给定的z；Step 9.1: Initialization: the number of iterations j=1,

According to the formula

Calculate the beamforming weight value

in

R_c (f_l ) is the covariance matrix of the clutter c[l]; in each iteration, a set of lattice point selection vectors {z₁ ,...,z_α } are randomly generated, for a given z ;步骤9.2：重复步骤9.3到步骤9.6的迭代过程：Step 9.2: Repeat the iterative process from Step 9.3 to Step 9.6:步骤9.3：随机产生一组格点选择矢量{z₁,...,z_α}；Step 9.3: Randomly generate a set of lattice point selection vectors {z₁ ,...,z_α };步骤9.4：根据公式w_g,l＝z⊙ξ_g,l计算

并构成

Step 9.4: Calculate according to the formula w_g,l =z⊙ξg_,l

and constitute

根据公式

计算r_g,l，将r_g,l表示成矢量r；According to the formula

Calculate r_g,l and represent r_g,l as a vector r;将

和r代入公式

得到x重构结果

和目标DOA估计结果

Will

and r into the formula

get x reconstruction result

and target DOA estimation results

根据公式

得到BMSE

According to the formula

get BMSE

步骤9.5：基于最小BMSE得到

Step 9.5: Obtain based on minimum BMSE

步骤9.6：根据

基于

更新

中相应的权重值得到

并令j＝j+1；其中

为Dc[l]的协方差矩阵，

其中，

z_op为格点选择矢量，c[l]为杂波；Step 9.6: According to

based on

renew

The corresponding weight values in the

and let j=j+1; where

is the covariance matrix of Dc[l],

in,

z_op is the grid selection vector, and c[l] is the clutter;步骤9.7：当

时，迭代停止，输出最优天线选择

e₀为事先设定的阈值。Step 9.7: When

When , the iteration stops and the optimal antenna selection is output

e₀ is a preset threshold.

Images