CN104076332A - 一种雷达均匀线性阵列幅度和相位的估计方法 - Google Patents

Abstract

本发明公开了一种雷达均匀线性阵列幅度和相位的估计方法，涉及阵列信号处理领域。包括：步骤1，设定雷达天线阵列具有M个阵元，其中前g个阵元的幅度和相位不存在误差；步骤2，将幅度、相位的预估计值作为均匀线性阵列的幅度、相位扰动初值；通过幅相扰动初值建立均匀线性阵列的阵列流型；步骤3，构建观测矩阵的子空间；步骤4，利用阵列流型及观测矩阵的子空间求解目标稀疏矩阵；步骤5，利用阵列流型及观测矩阵的子空间以及目标稀疏矩阵对阵列流型进行优化；步骤6，由优化后阵列流型矩阵对角线元素求出均匀线性阵列的幅度估计值和相位估计值。本发明在低信噪比且快拍数量较少的情况下，有效的对均匀线性阵列的幅度、相位进行估计。

Description

一种雷达均匀线性阵列幅度和相位的估计方法

技术领域

本发明属于雷达阵列信号处理领域，涉及一种雷达均匀线性阵列幅度和相位的估计方法。

背景技术

阵列信号处理是现代信号处理的一个重要研究分支。但在实际工程应用中，由于各种不可避免的误差，实际阵列流型往往会出现一定程度的偏差或扰动，这样，高分辨率的谱估计计算性能便会严重恶化，因而阵列的误差估计是阵列信号处理中很重要的研究方向之一。

由于各种扰动、误差的影响，当实际模型与假设模型不符合时，许多基于理想模型基础上的谱估计性能将严重下降，因而必须采取补偿或者阵列误差校准的措施。实际工程应用时，通常会出现快拍数量有限且信噪比低的问题。当快拍数少时或信噪比低时，信号与噪声没有足够的时间解相关，噪声协方差矩阵也没有收敛，因此会对后期各种参数估计性能造成影响。

目前，大多数阵列误差校准方法是通过对阵列扰动进行建模，将阵列误差校准逐渐转化为一个参数估计的问题的思想实现的。参数类的阵列校准方法可分为有源校准和自校准两类；有源校准通过在空间设置方位已知的辅助信源对阵列扰动参数进行离线估计，但这类校准算法对辅助信号源有较高的方位信息的要求，所以当辅助信号的方位信息有偏差时，这类校准算法会带来较大的偏差；而自校准通常根据某种优化函数对空间信源的方位与阵列的扰动参数进行联合估计，可以不需要方位已知的辅助信源，可以在线完成实际参数估计，所以校准的度较高，但由于误差参数与方位参数之间的耦合和某些病态的阵列结构，参数估计的唯一辨识往往无法保证，且其参数联合估计对应的高维、多模非线性优化问题带来了庞大的运算量。

在期刊IEEE Trans中，Friedlander B和Weiss A J基于子空间原理，提出一种信源方位、阵列互耦、阵元增益及相位扰动交替迭代估计的阵列自校准技术，但该技术需要求解高维非线性的优化问题，运算量庞大、收敛速度慢，全局收敛性无法保证，而且对于均匀线阵，阵列扰动的参数的估计存在模糊问题。在期刊IEEE Trans中，See C M S和Poh B K通过对阵列观测空间内一系列离散方位上阵列导向矢量的测量，基于最大似然准则，提出了一种信源方位、校准矩阵以及阵元位置同时联合优化估计的方法。虽然该方法通过参数联合估计和高斯—牛顿梯度算法与遗传算法相混合的优化算法来提高参数估计的收敛性和全局搜索性能，但其对应的非线性优化问题有很高的维数导致了更庞大的运算量。

《空间谱估计理论与算法》一书中提出了用辅助阵元法(ISM)来进行幅度相位误差的估计，通过引入少量校准的辅助阵元对阵元幅度相位误差进行无模糊联合估计，该方法只需要参数的一维搜索，不存在参数联合估计的局部收敛问题。但当快拍数据有限、信噪比低时，该估计结果会出现一些偏差。

发明内容

本发明的目的在于针对上述已有技术的不足，提出了一种雷达均匀线性阵列幅度和相位的估计方法，实现了在低信噪比且快拍数量较少的情况下，有效的对均匀线性阵列的幅度值、相位值进行估计，具有良好的稳健性。

为达到上述目的，本发明采用以下技术方案予以实现。

一种雷达均匀线性阵列幅度和相位的估计方法，其特征在于，包括以下步骤：

步骤1，设定雷达天线阵列为均匀线性阵列，均匀线性阵列具有M个阵元，其中前g个阵元的幅度和相位不存在误差；

步骤2，确定均匀线性阵列的导向矢量L和均匀线性阵列的幅度相位扰动矢量Γ；均匀线性阵列的导向矢量L表达式为：L＝[l₁l₂,...,l_h,...,l_M]，其中，l_h表示均匀线性阵列第h个阵元的导向矢量，1≤h≤M，均匀线性阵列的幅度相位扰动矢量Γ表达式Γ＝[1 1 1,...,Γ_M-g-1 Γ_M-g,...,Γ_M]^T，其中，扰动矢量Γ中前g个1代表幅度和相位不存在误差的阵元，后M-g个元素代表幅度和相位存在误差的阵元；M表示均匀线性阵列阵元数目，g表示幅度和相位不存在误差的阵元数目；

对均匀线性阵列的幅度相位扰动矢量Γ进行预估计，得到均匀线性阵列的幅度相位扰动预估计矢量Γ′；均匀线性阵列的幅度相位扰动预估计矢量Γ′表达式为：Γ′＝[1 1 1,...,Γ′_M-g-1 Γ′_M-g,...,Γ′_M]^T；(·)^T表示转置操作，扰动矢量Γ′中前g个1代表幅度和相位不存在误差的阵元，后M-g个元素代表幅度和相位存在误差的阵元；

利用预估计矢量Γ′和均匀线性阵列的导向矢量L构建均匀线性阵列的阵列流型D，即D＝diag(Γ′)·L，diag(.)表示矢量对角化，·表示矩阵点乘；

步骤3，雷达天线接收T个快拍的目标回波数据，利用目标回波数据以及阵列流型D建立观测矩阵Y^*；通过观测矩阵Y^*求解观测矩阵Y^*的子空间X；

步骤4，利用阵列流型D构造矩阵Φ＝DA(θ)，A(θ)是目标的导向矢量；再利用观测矩阵Y^*的子空间X求取矩阵Φ的支撑集Ω；构建支撑集Ω对应矩阵Φ中不为零的列向量组成的矩阵Φ_Ω，并利用观测矩阵Y^*的子空间X构建矩阵S_Ω；矩阵S_Ω组成目标稀疏矩阵S的非零列，目标的稀疏矩阵S的其余列为零，得到目标的稀疏矩阵S；

步骤5，利用观测矩阵Y^*的子空间X和目标的稀疏矩阵S对阵列流型D进行优化，得到优化后的阵列流型D^*；

步骤6，根据优化后的阵列流型D^*求取均匀线性阵列的幅度估计值和相位估计值，即对优化后的阵列流型D^*对角线上的数据求绝对值，得到均匀线性阵列的幅度估计值，对优化后的阵列流型D^*对角线上的数据求角度，即求取对角线上的数据的实部与虚部的反正切值，得到均匀线性阵列的相位估计值。

上述技术方案的特点和进一步改进在于：

(1)步骤3包括以下子步骤：

3a)接收T个快拍的目标回波数据，利用目标回波数据和阵列流型D建立观测矩阵Y^*：

Y^*＝DA(θ)S+E

其中，观测矩阵Y^*∈C^M×T,噪声矩阵E∈C^M×T，稀疏矩阵S∈C^N×T，A(θ)是目标的导向矢量，其中D表示阵列流型；M为该均匀线性阵列的阵元数，T为快拍数，N为采样点数，C为基矩阵；在[-90°,90°]间每隔1°进行采样，则采样点数N＝181；

3b)设定X＝U_SW^1/2为观测矩阵Y^*的子空间；

其中，W＝(Λ_S-σ_e ²I_K)²Λ_S ^-1为渐进最佳权值，U_S表示观测矩阵Y^*的K个目标的奇异值对应的奇异向量组成的矩阵，Λ_S为K个目标的奇异值组成的对角矩阵，且1≤K≤M；M表示均匀线性阵列阵元数目；λ_P表示观测矩阵Y^*的第P个噪声的奇异值，且1≤P≤M；σ_e ²表示噪声功率；I_K为K×K的单位矩阵，(·)^-1表示矩阵取逆操作。

(2)步骤4包括以下子步骤：

4a)设定矩阵Φ＝DA(θ)，A(θ)是目标的导向矢量，Ω为矩阵Φ的支撑集，支撑集Ω中包含了矩阵Φ中不为零的列的标号，元素β是支撑集Ω中的一个元素，通过求解以下公式得到第i次迭代中的元素β_i：

${\mathrm{\β}}_i=\arg \underset{p}{\min }\left[\underset{q}{\mathrm{\Σ}}\right|{\left({\mathrm{\Φ}}_{{\mathrm{\Ω}}_{i-1}}^H{R}_{i-1}\right|\left|{\mathrm{\Φ}}_{{\mathrm{\Ω}}_{i-1}}^H{R}_{i-1}{\left|\right|}_{\∞}\right)}_{\mathrm{pq}}\left|\right]$

其中，表示在第i-1次迭代中由支撑集Ω_i-1对应于矩阵Φ中不为零的列向量组成的矩阵，设定R_i-1表示在第i-1次迭代时观测矩阵Y^*的子空间X的矩阵，设定R₀为观测矩阵Y^*的子空间X；(·)^H表示共轭转置操作，p表示目标的稀疏矩阵S的行的序号，0<p≤N，N为采样点数；q表示目标的稀疏矩阵S的列的序号，0<q≤T，T为快拍数；arg(.)表示求解最优化，i＝1,2,3，…,K，K为目标数，||·||_∞表示无穷范数运算符；

4b)利用第i次迭代中的元素β_i求解第i次迭代时支撑集Ω_i＝Ω_i-1∪β_i；

4c)利用第i次迭代时支撑集Ω_i和观测矩阵Y^*的子空间X得到第i次迭代时矩阵R_i，表达式为下式：

$R_i=X-{\mathrm{\&Phi;}}_{{\mathrm{\&Omega;}}_i}{\left({\mathrm{\&Phi;}}_{{\mathrm{\&Omega;}}_i}^H{\mathrm{\&Phi;}}_{{\mathrm{\&Omega;}}_i}\right)}^{-1}{\mathrm{\&Phi;}}_{{\mathrm{\&Omega;}}_i}^{H}X,$

其中，表示在第i次迭代中由支撑集Ω_i对应于矩阵Φ中不为零的列向量组成的矩阵；X为观测矩阵Y^*的子空间；(·)^H表示共轭转置操作，(·)^-1表示矩阵取逆操作；

4d)令迭代次数i增加1，重复迭代以上几步4a)-4c)，直到i等于K，得到Ω_K；设定支撑集Ω＝Ω_K，求解出支撑集Ω；

4e)构建支撑集Ω对应于矩阵Φ中不为零的列向量组成的矩阵Φ_Ω，再利用观测矩阵Y^*的子空间X构建矩阵S_Ω，矩阵S_Ω公式为

其中，Φ_Ω表示支撑集Ω对应于矩阵Φ中不为零的列向量组成的矩阵，(·)^H表示共轭转置操作，X为观测矩阵Y^*的子空间，(·)^-1表示矩阵取逆操作；

4f)矩阵S_Ω组成目标稀疏矩阵S的非零列，目标的稀疏矩阵S的其余列为零，得到目标的稀疏矩阵S。

(3)步骤5包括以下子步骤：

5a)利用观测矩阵Y^*的子空间X，将均匀线性阵列的幅度相位估计转换为优化问题，也就是公式(1)：

$\underset{S,D}{\min }\left|\right|X-\mathrm{DA}\left(\mathrm{\&theta;}\right)S{\left|\right|}_F^2+\mathrm{\&lambda;}\left|\right|S{\left|\right|}_{\&infin;,0}---\left(1\right)$

其中，||·||_F为Frobenius范数运算符，X为观测矩阵Y^*的子空间，A(θ)是目标的导向矢量，S为目标的稀疏矩阵，D为均匀线性阵列的阵列流型，||S||_∞，0定义为||·||_∞，0表示混合范数，p表示目标的稀疏矩阵S的行的序号，0<p≤N，N为采样点数；q表示目标的稀疏矩阵S的列的序号，0<q≤T，T为快拍数，λ>0为正则化参数；

5b)将目标的稀疏矩阵S代入优化公式(1)得到优化求解阵列流型D的公式(2)；

$\underset{D}{\min }\left|\right|X-\mathrm{DA}\left(\mathrm{\&theta;}\right)S{\left|\right|}_F^2$

s.t.||S||_∞，0＝常数(约束条件)   (2)

5c)将优化求解阵列流型D的公式(2)进行矢量化，得到阵列流型D的矢量化形式d，对阵列流型D的矢量化形式d进行迭代，得到阵列流型D的矢量化形式d的第j次迭代值d_j的表达式为下式(3)：

$d_j=\arg \underset{d_{j-1}}{\min }\left|\right|\mathrm{vec}\left(X\right)-d_{j-1}A\left(\mathrm{\&theta;}\right)S{\left|\right|}_2^2---\left(3\right)$

其中，vec(·)表示矢量化，d₀为阵列流型D对角线元素，j表示迭代次数，||·||₂表示2范数运算符，X为观测矩阵Y^*的子空间，A(θ)是目标的导向矢量，S表示目标的稀疏矩阵，arg(.)表示求解最优化；

5d)将第j次迭代值d_j进行1范数加权得到第j次迭代值范数加权值d^* _j，d^* _j的表达式如下：

${d^{*}}_j=\arg \underset{{d^{*}}_{j-1}}{\min }\left|\right|\stackrel{\&OverBar;}{W}\mathrm{vec}\left(X\right)-{d^{*}}_{j-1}A\left(\mathrm{\&theta;}\right)S){\left|\right|}_1---\left(4\right)$

其中，W表示观测矩阵Y^*中行的互协方差矩阵的非零对角矩阵，其中I_M表示M×M的单位矩阵，M表示均匀线性阵列阵元数目，vec(·)表示矢量化，||·||₁表示1范数运算符，A(θ)是目标的导向矢量，S表示目标的稀疏矩阵，d^* ₀由d₀加权1范数得到，arg(.)表示求解最优化；

5e)求解公式(4)，若达到终止条件||d^* _j-d^* _j-1||₂<ε,取ε＝10^-4，设定第j次迭代值范数加权值d^* _j为阵列流型D矢量化并加权1范数矢量d^*；若没有达到终止条件，则令j增加1，继续求解公式(4)；

5f)通过阵列流型D矢量化并加权1范数矢量d^*得到优化后的阵列流型D^*，即D^*＝diag(d^*)；diag(.)表示矢量对角化。

与现有技术相比，本发明具有突出的实质性特点和显著的进步。

第一，本发明将信号子空间拟合和加权L1范数逼近的思想结合，发明的迭代算法可以对均匀线性阵列的幅度值和相位值进行估计。

第二，本发明的方法可以在快拍数量有限且低信噪比的情况下，降低均匀线性阵列幅度值和相位值估计误差。

本发明将子空间拟合和加权L1范数逼近结合的思想来进行幅度值和相位值的估计，主要应用于低信噪比、快拍数量有限的情况下均匀线性阵列幅度值和相位值的估计。

附图说明

下面结合附图和具体实施方式对本发明做进一步说明。

图1是本发明的实现流程图；

图2是本发明应用场景图；

图3是本发明迭代算法L1-WSSF估计幅度与仿真中设定的阵元真实幅度的对比图，横坐标为阵元数，纵坐标为幅度估计；

图4是本发明迭代算法L1-WSSF估计相位与仿真中设定的阵元真实相位的对比图，横坐标为阵元数，纵坐标为相位估计，单位度；

图5是本发明L1-WSSF法估计的阵列的幅度误差与现有技术ISM估计的幅度误差的对比图；横坐标为信噪比，单位dB；纵坐标为幅度估计误差；

图6是本发明L1-WSSF法估计的阵列的相位误差与现有技术ISM估计的相位误差的对比图；横坐标为信噪比，单位dB；纵坐标为相位估计误差。

具体实施方式

参照图1，说明本发明的一种雷达均匀线性阵列幅度和相位的估计方法，其具体实现步骤如下：

步骤1，设定雷达天线阵列为均匀线性阵列(Uniform Linear Array，ULA)且具有M个阵元，其中前g个阵元的幅度和相位不存在误差。

通过步骤1建立实验平台。

步骤2，确定均匀线性阵列的导向矢量L和均匀线性阵列的幅度相位扰动矢量Γ；对均匀线性阵列的幅度相位扰动矢量Γ进行预估计，得到均匀线性阵列的幅度相位扰动预估计矢量Γ′；利用预估计矢量Γ′和均匀线性阵列的导向矢量L构建均匀线性阵列的阵列流型D，即D＝diag(Γ′)·L，diag(.)表示矢量对角化，·表示矩阵点乘。

均匀线性阵列的导向矢量L表达式为：L＝[l ₁l₂,...,l_h,...,l_M]，其中，l_h表示均匀线性阵列第h个阵元的导向矢量，1≤h≤M，M为均匀线性阵列阵元数目；

均匀线性阵列的幅度相位扰动矢量Γ表达式为：

Γ＝[1 1 1,...,Γ_M-g-1 Γ_M-g,...,Γ_M]^T，其中，扰动矢量Γ中前g个1代表幅度和相位不存在误差的阵元，后M-g个元素代表幅度和相位存在误差的阵元；其中，M表示均匀线性阵列阵元数目，g表示幅度和相位不存在误差的阵元数目。

均匀线性阵列的幅度相位扰动预估计矢量Γ′表达式为：

Γ′＝[1 1 1,...,Γ′_M-g-1 Γ′_M-g,...,Γ′_M]^T；(·)^T表示转置操作；其中，扰动矢量Γ′中前g个1代表幅度和相位不存在误差的阵元，后M-g个元素代表幅度和相位存在误差的阵元。

在对均匀线性阵列幅度相位扰动矢量Γ进行预估计的过程中，使用辅助阵元法(instrumental sensor method，ISM)实现预估计，辅助阵元法出自王永良等编著的《空间谱估计理论与算法》一书，通过引入少量校准的阵元，利用信号子空间与噪声子空间正交的原理预估计出均匀线性阵列的幅度相位误差。

步骤3，雷达天线接收T个快拍的目标回波数据，利用目标回波数据以及阵列流型D建立观测矩阵Y^*；通过观测矩阵Y^*求解观测矩阵Y^*的子空间X。

3a)接收T个快拍的目标回波数据，利用目标回波数据和阵列流型D建立观测矩阵Y^*：

Y^*＝DA(θ)S+E

其中，观测矩阵Y^*∈C^M×T,噪声矩阵E∈C^M×T，稀疏矩阵S∈C^N×T，A(θ)是目标的导向矢量，其中D表示阵列流型；M为该均匀线性阵列的阵元数，T为快拍数，N为采样点数，C为基矩阵；在[-90°,90°]间每隔1°进行采样，则采样点数N＝181。

3b)设定X＝U_SW^1/2为观测矩阵Y^*的子空间；

其中，W＝(Λ_S-σ_e ²I_K)²Λ_S ^-1为渐进最佳权值，U_S表示观测矩阵Y^*的K个目标的奇异值对应的奇异向量组成的矩阵，Λ_S为K个目标的奇异值组成的对角矩阵，且1≤K≤M；M表示均匀线性阵列阵元数目；λ_P表示观测矩阵Y^*的第P个噪声的奇异值，且1≤P≤M；σ_e ²表示噪声功率；I_K为K×K的单位矩阵，(·)^-1表示矩阵取逆操作。

步骤4，利用阵列流型D构造矩阵Φ＝DA(θ)，A(θ)是目标的导向矢量；再利用观测矩阵Y^*的子空间X求取矩阵Φ的支撑集Ω；构建支撑集Ω对应于矩阵Φ中不为零的列向量组成的矩阵Φ_Ω，并利用观测矩阵Y^*的子空间X构建矩阵S_Ω；矩阵S_Ω组成目标稀疏矩阵S的非零列，目标的稀疏矩阵S的其余列为零，得到目标的稀疏矩阵S。

4a)设定矩阵Φ＝DA(θ)，A(θ)是目标的导向矢量，Ω为矩阵Φ的支撑集，支撑集Ω中包含了矩阵Φ中不为零的列的标号，元素β是支撑集Ω中的一个元素，通过求解以下公式得到第i次迭代中的元素β_i：

${\mathrm{\&beta;}}_i=\arg \underset{p}{\min }\left[\underset{q}{\mathrm{\&Sigma;}}\right|{\left({\mathrm{\&Phi;}}_{{\mathrm{\&Omega;}}_{i-1}}^H{R}_{i-1}\right|\left|{\mathrm{\&Phi;}}_{{\mathrm{\&Omega;}}_{i-1}}^H{R}_{i-1}{\left|\right|}_{\&infin;}\right)}_{\mathrm{pq}}\left|\right]$

其中，表示在第i-1次迭代中由支撑集Ω_i-1对应于矩阵Φ中不为零的列向量组成的矩阵，设定R_i-1表示在第i-1次迭代时观测矩阵Y^*的子空间X的矩阵，设定R₀为观测矩阵Y^*的子空间X；(·)^H表示共轭转置操作，p表示目标的稀疏矩阵S的行的序号，0<p≤N，N为采样点数；q表示目标的稀疏矩阵S的列的序号，0<q≤T，T为快拍数；arg(.)表示求解最优化，i＝1,2,3，…,K，K为目标数，||·||_∞表示无穷范数运算符。

4b)利用第i次迭代中的元素β_i求解第i次迭代时支撑集Ω_i＝Ω_i-1∪β_i；

4c)利用第i次迭代时支撑集Ω_i和观测矩阵Y^*的子空间X得到第i次迭代时矩阵R_i，表达式为下式：

$R_i=X-{\mathrm{\&Phi;}}_{{\mathrm{\&Omega;}}_i}{\left({\mathrm{\&Phi;}}_{{\mathrm{\&Omega;}}_i}^H{\mathrm{\&Phi;}}_{{\mathrm{\&Omega;}}_i}\right)}^{-1}{\mathrm{\&Phi;}}_{{\mathrm{\&Omega;}}_i}^{H}X,$

其中，表示在第i次迭代中由支撑集Ω_i对应于矩阵Φ中不为零的列向量组成的矩阵；X为观测矩阵Y^*的子空间；(·)^H表示共轭转置操作，(·)^-1表示矩阵取逆操作；

4d)令迭代次数i增加1，重复迭代以上几步4a)-4c)，直到i等于K，得到Ω_K；设定支撑集Ω＝Ω_K，求解出支撑集Ω；

4e)构建支撑集Ω对应于矩阵Φ中不为零的列向量组成的矩阵Φ_Ω，再利用观测矩阵Y^*的子空间X构建矩阵S_Ω，矩阵S_Ω公式为

其中，Φ_Ω表示支撑集Ω对应于矩阵Φ中不为零的列向量组成的矩阵，(·)^H表示共轭转置操作，X为观测矩阵Y^*的子空间，(·)^-1表示矩阵取逆操作；

4f)矩阵S_Ω组成目标稀疏矩阵S的非零列，目标的稀疏矩阵S的其余列为零，得到目标的稀疏矩阵S。

步骤5，利用观测矩阵Y^*的子空间X和目标的稀疏矩阵S对阵列流型D进行优化，得到优化后的阵列流型D^*。

5a)利用观测矩阵Y^*的子空间X，将均匀线性阵列的幅度相位估计转换为优化问题，也就是公式(1)：

$\underset{S,D}{\min }\left|\right|X-\mathrm{DA}\left(\mathrm{\&theta;}\right)S{\left|\right|}_F^2+\mathrm{\&lambda;}\left|\right|S{\left|\right|}_{\&infin;,0}---\left(1\right)$

其中，||·||_F为Frobenius范数运算符，X为观测矩阵Y^*的子空间，A(θ)是目标的导向矢量，S为目标的稀疏矩阵，D为均匀线性阵列的阵列流型，||S||_∞，0定义为||·||_∞，0表示混合范数，p表示目标的稀疏矩阵S的行的序号，0<p≤N，N为采样点数；q表示目标的稀疏矩阵S的列的序号，0<q≤T，T为快拍数，λ>0为正则化参数。

5b)将目标的稀疏矩阵S代入优化公式(1)得到优化求解阵列流型D的公式(2)；

$\underset{D}{\min }\left|\right|X-\mathrm{DA}\left(\mathrm{\&theta;}\right)S{\left|\right|}_F^2$

s.t.||S||_∞，0＝常数(约束条件)   (2)

现有技术的方法一般求解上式，则会遇到大量的数值错误，为了克服上述的缺点本发明中，不对上式直接求解，进行以下步骤：

5c)将优化求解阵列流型D的公式(2)进行矢量化，得到阵列流型D的矢量化形式d，对阵列流型D的矢量化形式d进行迭代，得到阵列流型D的矢量化形式d的第j次迭代值d_j的表达式为下式(3)：

$d_j=\arg \underset{d_{j-1}}{\min }\left|\right|\mathrm{vec}\left(X\right)-d_{j-1}A\left(\mathrm{\&theta;}\right)S{\left|\right|}_2^2---\left(3\right)$

其中，vec(·)表示矢量化，d₀为阵列流型D对角线元素，j表示迭代次数，||·||₂表示2范数运算符，X为观测矩阵Y^*的子空间，A(θ)是目标的导向矢量，S表示目标的稀疏矩阵，arg(.)表示求解最优化。

对于求解(3)式，当快拍数少，且噪声为有色噪声时，该表达式的稳健性就急剧下降。针对该问题，我们进行采取范数加权的方法解决。

5d)将第j次迭代值d_j进行1范数加权得到第j次迭代值范数加权值d^* _j，d^* _j的表达式如下：

${d^{*}}_j=\arg \underset{{d^{*}}_{j-1}}{\min }\left|\right|\stackrel{\&OverBar;}{W}\mathrm{vec}\left(X\right)-{d^{*}}_{j-1}A\left(\mathrm{\&theta;}\right)S){\left|\right|}_1---\left(4\right)$

其中，W表示观测矩阵Y^*中行的互协方差矩阵的非零对角矩阵，其中I_M表示M×M的单位矩阵，M表示均匀线性阵列阵元数目，vec(·)表示矢量化，||·||₁表示1范数运算符，A(θ)是目标的导向矢量，S表示目标的稀疏矩阵，d^* ₀由d₀加权1范数得到，arg(.)表示求解最优化。

在本发明中，当W≠I_M时，优化求解阵列流型D的方法称为子空间拟合和加权L1范数法—L1-WSSF(L1-norm weighted sparse subspace fitting)；当W＝I_M时，优化求解阵列流型D的方法称为子空间拟合和L1范数—L1-SSF(L1-norm sparse subspace fitting)。

5e)求解公式(4)，若达到终止条件||d^* _j-d^* _j-1||₂<ε,取ε＝10^-4，设定第j次迭代值范数加权值d^* _j为阵列流型D矢量化并加权1范数矢量d^*；若没有达到终止条件，则令j增加1，继续求解公式(4)。

5f)通过阵列流型D矢量化并加权1范数矢量d^*得到优化后的阵列流型D^*，即D^*＝diag(d^*)；diag(.)表示矢量对角化。

步骤6，根据优化后的阵列流型D^*求取均匀线性阵列的幅度估计值和相位估计值，即对优化后的阵列流型D^*对角线上的数据求绝对值，得到均匀线性阵列的幅度估计值，对优化后的阵列流型D^*对角线上的数据求角度，即求取对角线上的数据的实部与虚部的反正切值，得到均匀线性阵列的相位估计值。

本发明的效果可以通过以下仿真结果进一步说明。

1.仿真条件：

本发明的仿真是在MATLAB R2009a的软件环境下进行的。

2.仿真内容：

(1)L1-WSSF法得到的均匀线性阵列的幅度估计值和相位估计值分析。

仿真环境使用图2所示的雷达阵地图，采用3根幅度和相位不存在误差的辅助天线和10根幅度和相位存在误差的天线构成均匀线性阵列(ULA)。假设有空间存在两个目标，到达均匀线性阵列的波达方向分别为θ₁＝10°和θ₂＝-20°，快拍数为10个，信噪比SNR为50dB,对环境中的角度[-90°,90°]按照间隔1°进行划分。如图3所示，将本发明迭代算法L1-WSSF得到的均匀线性阵列的幅度估计值与仿真中设定的均匀线性阵列的幅度值进行对比，横坐标为阵元数，纵坐标为幅度估计；其中，图中所示的真实幅度即为仿真中设定的均匀线性阵列的幅度值；如图4所示，将本发明迭代算法L1-WSSF得到的均匀线性阵列的相位估计值与仿真中设定的的均匀线性阵列的相位值进行对比，横坐标为阵元数，纵坐标为相位估计；其中，图中所示的真实相位即为仿真中设定的均匀线性阵列的相位值。

从图3、图4可看出本发明在快拍数量(10个快拍)有限的条件下，能正确得到均匀线性阵列的幅度估计值与相位估计值。现有技术中的辅助阵元法ISM求解均匀线性阵列幅度估计值和相位估计值时(也就是对步骤2得出阵列流型D对角线上的数据求绝对值，得到ISM法求解出的均匀线性阵列的幅度估计值，对阵列流型D对角线上的数据求角度，即求取对角线上的数据的实部与虚部的反正切值，得到ISM法求解出的均匀线性阵列的相位估计值)，通常需要使用几十上百个快拍数。

(2)本发明L1-WSSF法与现有技术辅助阵元法ISM法求解均匀线性阵列的幅度估计值和相位估计值性能的比较。

采用100次Monto Carlo实验，将稀疏子空间拟合和加权L1范数法(L1-WSSF)与辅助阵元法(ISM)做比较，其中定义幅度误差为 表示均匀线性阵列的幅度估计值，a₁,a₂…a₁₃为仿真中设定的均匀线性阵列的幅度值；定义相位误差为表示均匀线性阵列相位估计值，为仿真中设定的均匀线性阵列的相位值，快拍数为10，仿真结果如图5、图6所示。

仿真条件中快拍数仍为10个，从图5横坐标为信噪比，纵坐标为幅度估计误差，从图5可看出本发明L1-WSSF法得到的阵列幅度估计误差小于现有技术ISM得到的幅度估计误差；图6横坐标为信噪比，纵坐标为相位估计误差。从图6可看出本发明L1-WSSF法得到的阵列相位估计误差小于现有技术ISM得到的相位估计误差；可以看出本发明方法可以较准确、稳健的得到均匀线性阵列的幅度估计值和相位估计值。在快拍数有限、低信噪比的情况下，较之辅助阵元法(ISM)均方误差更低，性能更好。

Claims (4)

1.一种雷达均匀线性阵列幅度和相位的估计方法，其特征在于，包括以下步骤：

步骤1，设定雷达天线阵列为均匀线性阵列，均匀线性阵列具有M个阵元，其中前g个阵元的幅度和相位不存在误差；

步骤2，确定均匀线性阵列的导向矢量L和均匀线性阵列的幅度相位扰动矢量Γ；均匀线性阵列的导向矢量L表达式为：L＝[l₁l₂,...,l_h,...,l_M]，其中，l_h表示均匀线性阵列第h个阵元的导向矢量，1≤h≤M，均匀线性阵列的幅度相位扰动矢量Γ表达式Γ＝[1 1 1,...,Γ_M-g-1 Γ_M-g,...,Γ_M]^T，其中，扰动矢量Γ中前g个1代表幅度和相位不存在误差的阵元，后M-g个元素代表幅度和相位存在误差的阵元；M表示均匀线性阵列阵元数目，g表示幅度和相位不存在误差的阵元数目；

对均匀线性阵列的幅度相位扰动矢量Γ进行预估计，得到均匀线性阵列的幅度相位扰动预估计矢量Γ′；均匀线性阵列的幅度相位扰动预估计矢量Γ′表达式为：Γ′＝[1 1 1,...,Γ′_M-g-1 Γ′_M-g,...,Γ′_M]^T；(·)^T表示转置操作，扰动矢量Γ′中前g个1代表幅度和相位不存在误差的阵元，后M-g个元素代表幅度和相位存在误差的阵元；

利用预估计矢量Γ′和均匀线性阵列的导向矢量L构建均匀线性阵列的阵列流型D，即D＝diag(Γ′)·L，diag(.)表示矢量对角化，·表示矩阵点乘；

步骤3，雷达天线接收T个快拍的目标回波数据，利用目标回波数据以及阵列流型D建立观测矩阵Y^*；通过观测矩阵Y^*求解观测矩阵Y^*的子空间X；

步骤4，利用阵列流型D构造矩阵Φ＝DA(θ)，A(θ)是目标的导向矢量；再利用观测矩阵Y^*的子空间X求取矩阵Φ的支撑集Ω；构建支撑集Ω对应于矩阵Φ中不为零的列向量的矩阵Φ_Ω，并利用观测矩阵Y^*的子空间X构建矩阵S_Ω；矩阵S_Ω组成目标稀疏矩阵S的非零列，目标的稀疏矩阵S的其余列为零，得到目标的稀疏矩阵S；

步骤5，利用观测矩阵Y^*的子空间X和目标的稀疏矩阵S对阵列流型D进行优化，得到优化后的阵列流型D^*；

步骤6，根据优化后的阵列流型D^*求取均匀线性阵列的幅度估计值和相位估计值，即对优化后的阵列流型D^*对角线上的数据求绝对值，得到均匀线性阵列的幅度估计值，对优化后的阵列流型D^*对角线上的数据求角度，即求取对角线上的数据的实部与虚部的反正切值，得到均匀线性阵列的相位估计值。

2.根据权利要求1所述的一种雷达均匀线性阵列幅度和相位的估计方法，其特征在于，步骤3包括以下子步骤：

3a)接收T个快拍的目标回波数据，利用目标回波数据和阵列流型D建立观测矩阵Y^*：

Y^*＝DA(θ)S+E

其中，观测矩阵Y^*∈C^M×T,噪声矩阵E∈C^M×T，稀疏矩阵S∈C^N×T，A(θ)是目标的导向矢量，其中D表示阵列流型；M为该均匀线性阵列的阵元数，T为快拍数，N为采样点数，C为基矩阵；在[-90°,90°]间每隔1°进行采样，则采样点数N＝181；

3b)设定X＝U_SW^1/2为观测矩阵Y^*的子空间；

其中，W＝(Λ_S-σ_e ²I_K)²Λ_S ^-1为渐进最佳权值，U_S表示观测矩阵Y^*的K个目标的奇异值对应的奇异向量组成的矩阵，Λ_S为K个目标的奇异值组成的对角矩阵，且1≤K≤M；M表示均匀线性阵列阵元数目；λ_P表示观测矩阵Y^*的第P个噪声的奇异值，且1≤P≤M；σ_e ²表示噪声功率；I_K为K×K的单位矩阵，(·)^-1表示矩阵取逆操作。

3.根据权利要求1所述的一种雷达均匀线性阵列幅度和相位的估计方法，其特征在于，步骤4包括以下子步骤：

4a)设定矩阵Φ＝DA(θ)，A(θ)是目标的导向矢量，Ω为矩阵Φ的支撑集，支撑集Ω中包含了矩阵Φ中不为零的列的标号，元素β是支撑集Ω中的一个元素，通过求解以下公式得到第i次迭代中的元素β_i：

${\mathrm{\&beta;}}_i=\arg \underset{p}{\min }\left[\underset{q}{\mathrm{\&Sigma;}|}{\left({\mathrm{\&Phi;}}_{{\mathrm{\&Omega;}}_{i-1}}^H{R}_{i-1}\right|\left|{\mathrm{\&Phi;}}_{{\mathrm{\&Omega;}}_{i-1}}^H{R}_{i-1}{\left|\right|}_{\&infin;}\right)}_{\mathrm{pq}}\right|]$

其中，表示在第i-1次迭代中由支撑集Ω_i-1对应于矩阵Φ中不为零的列向量组成的矩阵，设定R_i-1表示在第i-1次迭代时观测矩阵Y^*的子空间X的矩阵，设定R₀为观测矩阵Y^*的子空间X；(·)^H表示共轭转置操作，p表示目标的稀疏矩阵S的行的序号，0<p≤N，N为采样点数；q表示目标的稀疏矩阵S的列的序号，0<q≤T，T为快拍数；arg(.)表示求解最优化，i＝1,2,3，…,K，K为目标数，||·||_∞表示无穷范数运算符；

4b)利用第i次迭代中的元素β_i求解第i次迭代时支撑集Ω_i＝Ω_i-1∪β_i；

4c)利用第i次迭代时支撑集Ω_i和观测矩阵Y^*的子空间X得到第i次迭代时矩阵R_i表达式为下式：

$R_i=X-{\mathrm{\&Phi;}}_{{\mathrm{\&Omega;}}_i}{\left({\mathrm{\&Phi;}}_{{\mathrm{\&Omega;}}_i}^H{\mathrm{\&Phi;}}_{{\mathrm{\&Omega;}}_i}\right)}^{-1}{\mathrm{\&Phi;}}_{{\mathrm{\&Omega;}}_i}^{H}X,$

其中，表示在第i次迭代中由支撑集Ω_i对应于矩阵Φ中不为零的列向量组成的矩阵；X为观测矩阵Y^*的子空间；(·)^H表示共轭转置操作，(·)^-1表示矩阵取逆操作；

4d)令迭代次数i增加1，重复迭代以上几步4a)-4c)，直到i等于K，得到Ω_K；设定支撑集Ω＝Ω_K，求解出支撑集Ω；

4e)构建支撑集Ω对应于集合Φ中不为零的列向量的矩阵Φ_Ω，再利用观测矩阵Y^*的子空间X构建矩阵S_Ω，矩阵S_Ω公式为

其中，Φ_Ω表示支撑集Ω对应于矩阵Φ中的列向量集合，(·)^H表示共轭转置操作，X为观测矩阵Y^*的子空间，(·)^-1表示矩阵取逆操作；

4f)矩阵S_Ω组成目标稀疏矩阵S的非零列，目标的稀疏矩阵S的其余列为零，得到目标的稀疏矩阵S。

4.根据权利要求1所述的一种雷达均匀线性阵列幅度和相位的估计方法，其特征在于，步骤5包括以下子步骤：

5a)利用观测矩阵Y^*的子空间X，将均匀线性阵列的幅度相位估计转换为优化问题，也就是公式(1)：

$\underset{S,D}{\min }\left|\right|X-\mathrm{DA}\left(\mathrm{\&theta;}\right)S{\left|\right|}_F^2+\mathrm{\&lambda;}\left|\right|S{\left|\right|}_{\&infin;,0}---\left(1\right)$

其中，||·||_F为Frobenius范数运算符，X为观测矩阵Y^*的子空间，A(θ)是目标的导向矢量，S为目标的稀疏矩阵，D为均匀线性阵列的阵列流型，||S||_∞，0定义为||·||_∞，0表示混合范数，p表示目标的稀疏矩阵S的行的序号，0<p≤N，N为采样点数；q表示目标的稀疏矩阵S的列的序号，0<q≤T，T为快拍数，λ>0为正则化参数；

5b)将目标的稀疏矩阵S代入优化公式(1)得到优化求解阵列流型D的公式(2)；

$\underset{D}{\min }\left|\right|X-\mathrm{DA}\left(\mathrm{\&theta;}\right)S{\left|\right|}_F^2$

s.t.||S||_∞，0＝常数(约束条件)   (2)

5c)将优化求解阵列流型D的公式(2)进行矢量化，得到阵列流型D的矢量化形式d，对阵列流型D的矢量化形式d进行迭代，得到阵列流型D的矢量化形式d的第j次迭代值d_j的表达式为下式(3)：

$d_j=\arg \underset{d_{j-1}}{\min }\left|\right|\mathrm{vec}\left(X\right)-d_{j-1}A\left(\mathrm{\&theta;}\right)S{\left|\right|}_2^2---\left(3\right)$

其中，vec(·)表示矢量化，d₀为阵列流型D对角线元素，j表示迭代次数，||·||₂表示2范数运算符，X为观测矩阵Y^*的子空间，A(θ)是目标的导向矢量，S表示目标的稀疏矩阵，arg(.)表示求解最优化；

5d)将第j次迭代值d_j进行1范数加权得到第j次迭代值范数加权值d^* _j，d^* _j的表达式如下：

${d^{*}}_j=\arg \underset{{d^{*}}_{j-1}}{\min }\left|\right|\stackrel{\&OverBar;}{W}\mathrm{vec}\left(X\right)-{d^{*}}_{j-1}A\left(\mathrm{\&theta;}\right)S){\left|\right|}_1---\left(4\right)$

其中，W表示观测矩阵Y^*中行的互协方差矩阵的非零对角矩阵，其中I_M表示M×M的单位矩阵，M表示均匀线性阵列阵元数目，vec(·)表示矢量化，||·||₁表示1范数运算符，A(θ)是目标的导向矢量，S表示目标的稀疏矩阵，d^* ₀由d₀加权1范数得到，arg(.)表示求解最优化；

5e)求解公式(4)，若达到终止条件||d^* _j-d^* _j-1||₂<ε,取ε＝10^-4，设定第j次迭代值范数加权值d^* _j为阵列流型D矢量化并加权1范数矢量d^*；若没有达到终止条件，则令j增加1，继续求解公式(4)；

5f)通过阵列流型D矢量化并加权1范数矢量d^*得到优化后的阵列流型D^*，即D^*＝diag(d^*)；diag(.)表示矢量对角化。