CN103885042A - 基于杂波子空间的阵元误差估计方法 - Google Patents

Abstract

本发明公开了一种基于杂波子空间的阵元误差估计方法，其实现过程为：1.利用雷达先验信息计算杂波空时导向矢量矩阵；2.分别对雷达接收数据和杂波空时导向矢量矩阵进行奇异值分解；3.计算主奇异值的数目；4.根据得到的主奇异值的数目，构造杂波正交补空间矩阵；5.利用得到的杂波正交补空间和对接收数据分解得到的最大奇异值对应的左奇异向量的正交性，估计天线阵元误差；6.估计出的天线阵元误差修正杂波信号模型，进行杂波抑制。本发明实现了在样本数较少和脉冲数较少的情况下，均能取得较好的天线阵元误差估计性能，使得本发明具有误差估计精度高和估计性能稳健的优点。

Description

基于杂波子空间的阵元误差估计方法

技术领域

本发明属于通信技术领域，更进一步涉及机载雷达技术领域中的一种基于杂波子空间的阵元误差估计方法。本发明利用估计出来的阵元误差对理论杂波模型进行修正，可以显著提高机载雷达的杂波抑制性能和目标检测性能。

背景技术

机载雷达在检测地面低速运动目标时，通常面临强杂波的竞争，给运动目标检测带来了严峻的挑战。空时自适应处理（space time adaptive processing，STAP）可以有效抑制杂波并实现运动目标检测。通常STAP中的协方差矩阵是通过选择待检测单元邻近的单元作为训练样本来估计得到的。在训练样本独立同分布且和待检测单元的统计特性一致的条件下，样本数目大于两倍系统自由度时自适应处理的性能损失相对于最优处理而言不超过3dB。但是实际的杂波场景往往呈现非均匀性，破坏了样本的独立同分布假设，从而影响了机载雷达的杂波抑制和目标检测性能。特别在地面运动目标检测（ground moving target indicator，GMTI）时，由于道路密集，运动目标较多，训练样本中常常包含有目标信号。

中国科学院电子学研究所拥有的专利技术“一种基于样本训练的机载雷达非平稳杂波抑制方法”（专利申请号：2011102289306.X，申请公告号：CN103018727A）中提出了一种基于样本训练的机载雷达非平稳杂波抑制方法。该专利技术在利用训练样本进行空时自适应处理（STAP）时，通过减少距离向的训练样本，和以滑动窗形式依次滑动选取脉冲，来增加脉冲维的训练样本，有效提高了杂波协方差矩阵的估计精度，改善主瓣区杂波抑制性能。但是该专利技术仍然存在的不足是，当训练样本中包含有目标信号时，利用这些样本构造协方差矩阵并计算自适应权值，将导致杂波的抑制性能和目标的检测性能变差。

中国科学院电子学研究所拥有的专利技术“基于噪声失效奇异值分解的地基雷达杂波抑制方法”（专利申请号：201110028151.4，申请公告号：CN102621534A）中提出了一种基于噪声失效奇异值分解的地基雷达杂波抑制方法。该专利技术通过选择所有门限值中使得噪声向量最接近高斯向量的那个值，作为噪声系数的门限值，该方法当噪声向量满足高斯分布时，有效改善了杂波抑制效果。但是该专利技术仍然存在的不足是，当噪声向量不满足高斯分布时，利用该方法进行杂波抑制，将导致杂波的抑制性能和目标的检测性能变差。

发明内容

本发明的目的在于克服上述现有技术的不足，提出一种基于杂波子空间的阵元误差估计方法。本发明利用雷达先验信息计算杂波空时导向矢量矩阵，再分别对接收数据和杂波空时导向矢量矩阵进行奇异值分解，分别得到最大奇异值对应的左奇异向量和杂波的正交补空间，最后利用两者的正交性，来估计天线阵元误差。

本发明的具体实施步骤如下：

(1)计算杂波空时导向矢量矩阵：

利用雷达先验信息，计算杂波空时导向矢量矩阵。

(2)奇异值分解：

2a)对雷达接收数据进行奇异值分解，得到接收数据左奇异值矢量矩阵，主、次奇异值对应的对角矩阵和接收数据右奇异值矢量矩阵；

2b)对杂波空时导向矢量矩阵进行奇异值分解，得到左奇异值矢量矩阵，主、次奇异值对应的对角矩阵和右奇异值矢量矩阵。

(3)计算主奇异值的数目：

按照下式，计算得到主奇异值的数目：

$r=\min (\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{\mathrm{\σ}}_i\≥0.999\underset{i=1}{\overset{\mathrm{NM}}{\mathrm{\Σ}}}{\mathrm{\σ}}_i)$

其中，r表示主奇异值的数目，m表示奇异值的数目，m的取值范围为：m＝1,…,NM，N表示天线阵元数，M表示雷达一个相干处理间隔内发射的脉冲数，i表示奇异值个数的标号，i的取值范围为：i＝1,…,m，σ_i表示对杂波空时导向矢量矩阵进行奇异值分解得到的奇异值，min(·)表示取最小值操作。

(4)构造正交补空间矩阵：

取出左奇异矢量矩阵，去掉左奇异矢量矩阵中主奇异值对应矩阵的列，将剩余列构成杂波的正交补空间矩阵。

(5)估计天线阵元误差：

5a)取出接收数据左奇异矢量矩阵中的首列，利用其与杂波的正交补空间矩阵的正交性，按照下式，估计天线阵元锥削矢量：

$\underset{\hat{t}}{\min }\left({\left|\right|U^H(e_1\·/\hat{t})\left|\right|}_2^2\right)$

其中，

表示估计的天线阵元锥削矢量，min(·)表示取最小值操作，U表示杂波的正交补空间矩阵，H表示矩阵共轭转置操作，e₁表示接收数据左奇异矢量矩阵中的首列，||·||₂表示2范数操作，·/表示对向量对应元素相除的点除操作；

5b)按照下式，估计天线阵元误差：

$\mathrm{\Δ}\hat{t}=\hat{t}(1:N)-1_N$

其中，

表示估计的天线阵元误差，

表示估计的天线阵元锥削矢量，N表示天线阵元数，1_N表示N个全为1的列矢量。

(6)杂波抑制：

6a)利用M个全为1的列矢量构成的理想时域导向矢量，与估计出的天线阵元误差进行Kronecker积运算，得到修正后的空时锥削矢量，利用此空时锥削矢量与杂波空时导向矢量进行Hadamard积运算，得到修正后的杂波信号模型；

6b)利用杂波信号模型建立杂波协方差矩阵，进行杂波抑制。

本发明与现有技术相比具有如下优点：

第一，由于本发明利用雷达先验知识计算得到杂波对应的正交补空间，再通过对接收数据分解，得到最大奇异值对应的左奇异向量，最后利用两者的正交性，估计天线阵元误差，克服了现有技术中样本数据中含有目标，导致杂波抑制效果变差的问题，使得本发明可以在样本中含有目标时，仍可实现较高的天线阵元误差估计精度。

第二，由于本发明利用正交性来估计天线阵元误差，克服了现有技术中需要较高的多普勒分辨率和较多的距离单元平均的问题，使得本发明在样本数较少和脉冲数较少的情况下，均能取得较好的误差估计性能。

第三，由于本发明采用数学不等式确定杂波正交补空间的维数，克服了现有技术中，需对估计模型做出假设的缺点，使得本发明具有精度高和误差小的优点。

附图说明

图1是本发明的流程图；

图2是本发明与主杂波方法随脉冲数变化时天线阵元误差估计的准确度示意图；

图3是本发明与主杂波方法随样本数变化时天线阵元误差估计的准确度示意图；

图4是天线阵元存在误差时本发明与广义内积方法STAP输出信杂噪比结果示意图。

具体实施方式

下面结合附图对本发明做进一步的描述。

参照图1，本发明的具体实施步骤如下：

步骤1.计算杂波空时导向矢量矩阵。

利用雷达先验信息，计算杂波空时导向矢量矩阵。

第一步，按照下式，计算归一化多普勒频率：

其中，w₁表示归一化多普勒频率，v表示雷达载机速度，θ表示地面杂波散射体的俯仰角，表示地面杂波散射体的方位角，λ表示雷达工作波长，f_r表示雷达脉冲发射的重复频率。

第二步，按照下式，计算归一化空间频率：

其中，w₂表示归一化空间频率，d表示天线阵元间距，θ表示地面杂波散射体的俯仰角，

表示地面杂波散射体的方位角，ψ表示载机速度方向与雷达天线阵列轴向的夹角，λ表示雷达工作波长。

第三步，按照下式，计算时域导向矢量：

$u\left(w_1\right)={\left[\begin{array}{ccccccc}1& e^{j2\mathrm{\π}{w}_1*1}& e^{j2\mathrm{\π}{w}_1*2}& \·\·\·& e^{j2\mathrm{\π}{w}_1*(m_1-1)}& \·\·\·& e^{j2\mathrm{\π}{w}_1(M-1)}\end{array}\right]}^T$

其中，u(w₁)表示时域导向矢量，w₁表示归一化多普勒频率，m₁表示归一化多普勒频率的数目，m₁的取值范围为：m₁＝1,…,M，M表示雷达一个相干处理间隔内发射的脉冲数，

表示第m₁个归一化多普勒频率对应的相位，T表示转置操作。

第四步，按照下式，计算空域导向矢量：

$u\left(w_2\right)={\left[\begin{array}{ccccccc}1& e^{j2\mathrm{\π}{w}_2*1}& e^{j2\mathrm{\π}{w}_2*2}& \·\·\·& e^{j2\mathrm{\π}{w}_2*(m_2-1)}& \·\·\·& e^{j2\mathrm{\π}{w}_2*(N-1)}\end{array}\right]}^T$

其中，u(w₂)表示空域导向矢量，w₂表示归一化空间频率，m₂表示归一化多普勒频率的数目，m₂的取值范围为：m₂＝1,…,N，表示第m₂个归一化空间频率对应的相位，N表示雷达天线阵元数。

第五步，按照下式，计算空时导向矢量矩阵：

$V=u\left(w_1\right)\⊗u\left(w_2\right)$

其中，V表示空时导向矢量矩阵，u(w₁)表示时域导向矢量，w₁表示归一化多普勒频率，u(w₂)表示空域导向矢量，w₂表示归一化空间频率，

表示Kronecker积。

步骤2.奇异值分解。

对雷达接收数据进行奇异值分解，得到接收数据左奇异值矢量矩阵，主、次奇异值对应的对角矩阵和接收数据右奇异值矢量矩阵。

对杂波空时导向矢量矩阵进行奇异值分解，得到左奇异值矢量矩阵，主、次奇异值对应的对角矩阵和右奇异值矢量矩阵。

步骤3.计算主奇异值的数目。

按照下式，计算得到主奇异值的数目：

$r=\min (\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{\mathrm{\σ}}_i\≥0.999\underset{i=1}{\overset{\mathrm{NM}}{\mathrm{\Σ}}}{\mathrm{\σ}}_i)$

其中，r表示主奇异值的数目，m表示奇异值的数目，m的取值范围为：m＝1,…,NM，N表示天线阵元数，M表示雷达一个相干处理间隔内发射的脉冲数，i表示奇异值个数的标号，i的取值范围为：i＝1,…,m，σ_i表示对杂波空时导向矢量矩阵进行奇异值分解得到的奇异值，min(·)表示取最小值操作。

步骤4.构造正交补空间矩阵。

取出左奇异矢量矩阵，去掉左奇异矢量矩阵中主奇异值对应矩阵的列，将剩余列构成杂波的正交补空间矩阵。

步骤5.估计天线阵元误差。

取出接收数据左奇异矢量矩阵中的首列，利用其与杂波的正交补空间矩阵的正交性，按照下式，估计天线阵元锥削矢量：

$\underset{\hat{t}}{\min }\left({\left|\right|U^H(e_1\·/\hat{t})\left|\right|}_2^2\right)$

其中，表示估计的天线阵元锥削矢量，min(·)表示取最小值操作，U表示杂波的正交补空间矩阵，H表示矩阵共轭转置操作，e₁表示接收数据左奇异矢量矩阵中的首列，||·||₂表示2范数操作，·/表示对向量对应元素相除的点除操作。

按照下式，估计天线阵元误差：

$\mathrm{\Δ}\hat{t}=\hat{t}(1:N)-1_N$

其中，

表示估计的天线阵元误差，

表示估计的天线阵元锥削矢量，N表示天线阵元数，1_N表示N个全为1的列矢量。

步骤6.杂波抑制。

利用M个全为1的列矢量构成的理想时域导向矢量，与估计出的天线阵元误差进行Kronecker积运算，得到修正后的空时锥削矢量，利用此空时锥削矢量与杂波空时导向矢量进行Hadamard积运算，得到修正后的杂波信号模型。

利用杂波信号模型建立杂波协方差矩阵，进行杂波抑制。

下面结合附图2、附图3、附图4对本发明的效果做进一步的描述：

附图2是以本发明与主杂波方法估计出来的天线阵元误差与实际天线阵元误差的夹角余弦作为衡量标准，对比本发明与主杂波方法估计出来的天线阵元误差随脉冲数变化时，天线阵元误差估计准确程度的示意图。图2天线阵元误差估计准确程度的仿真条件和仿真结果如下：

1.仿真条件：

运行系统为Intel(R)Core(TM)2Duo CPU E45002.20GHz，64位Windows操作系统，仿真软件采用MATLAB R（2010a），仿真参数设置如下表所示：

参数	数值
		载频/GHz	1.2
波长/m	0.25
		阵元间距/m	0.125
阵元数	8
		脉冲重复频率/Hz	3000
发射带宽/MHz	2
		载机高度/km	5
脉冲数	20、40、60、80、100
		杂噪比	40
采样率/MHz	2
		载机速度/m/s	150
主波束俯仰角/。	-3
		主波束方位角/。	0
载机的偏航角/。	0

仿真参数如上表所示，待检测单元选200号距离门。导向矢量误差的方差

设为0.01。为了验证样本数据含有目标的情况下，本发明的有效性，试验中加入了5个目标，信杂比SNR分别为20、15、10、10、10，目标分别位于样本点的1、30、50、70、90处。目标的方位角为0°，俯仰角为-3°，归一化多普勒频率为0.25。计算得到的构型参数主波束的俯仰角为0°，载机的速度为145m/s。

2.仿真结果：

图2表示本发明与主杂波方法随脉冲数变化时，天线阵元误差估计的准确程度的仿真结果示意图，图2横坐标表示脉冲数，纵坐标表示夹角余弦。图2中的两条曲线分别表示本发明和主杂波方法的天线阵元误差估计性能随脉冲数的变化关系。

由图2中两种方法形成的夹角余弦曲线可见，在样本数不变的情况下，脉冲数较少时，本发明和实际误差比较接近，而主杂波方法相差较大。随着脉冲数逐渐增多，主杂波方法对天线阵元误差估计的准确程度逐渐逼近于本发明。主杂波方法在低脉冲的情况下阵元误差估计性能较差是因为在脉冲数目较少的情况下，多普勒分辨率较差；随着脉冲数目的增多，多普勒分辨率提高，主瓣对应的多普勒带宽变窄，其性能得到改善。

附图3是以本发明与主杂波方法估计出来的天线阵元误差与实际天线阵元误差的夹角余弦作为衡量标准，对比本发明与主杂波方法估计出来的天线阵元误差随样本数变化时，天线阵元误差估计准确程度示意图。图3天线阵元误差估计准确程度仿真条件和仿真结果如下：

1.仿真条件：

运行系统和仿真软件同上，除脉冲数设定为80，恒定不变，样本数分别为20、40、60、80、100外，其他参数设置同上。

2.仿真结果：

图3表示本发明与主杂波方法随样本数变化时，天线阵元误差估计的准确程度的仿真结果图，图2横坐标表示样本数，纵坐标表示夹角余弦。图中的两条曲线分别表示本发明和主杂波方法的天线阵元误差估计性能随样本数的变化关系。

由图3中两种方法形成的夹角余弦曲线可见，在脉冲数不变的情况下，样本数较少时，本发明和实际误差比较接近，而主杂波方法相差较大。随着样本数逐渐增多，主杂波方法对天线阵元误差估计的准确程度逐渐逼近于本发明。主杂波方法在低样本的情况下阵元误差估计性能较差是因为在低样本的情况下，噪声对估计性能的影响较大，随着样本数目的增加，噪声起伏得到平均。而本发明在低样本和低脉冲的情况下均可以取得较好的性能，这是因为本发明假定杂波由一个低秩的子空间张成，利用先验知识可以精确计算出其对应的正交补空间，然后通过对接收数据进行奇异值分解，取出其中杂噪比最大的向量，再利用子空间正交这一特性，可以准确估计出天线阵元误差。

附图4是天线阵元存在误差时，本发明与基于知识广义内积方法STAP输出信杂噪比结果图，仿真条件和仿真结果如下：

1.仿真条件：

运行系统和仿真软件同上，仿真主要参数设置如下：

主波束俯仰角是-3度，主波束方位角是0度，待检测单元位于第100号距离单元，采用全空时自适应处理算法，训练样本选取待检测单元左右各1.5倍自由度的样本。在样本中加入20个目标信号，强弱目标各占一半，信噪比分别为15dB与5dB，目标角度位于天线主瓣内，目标归一化多普勒频率为0.2。

仿真基本参数设置如下表所示：

参数	数值
		载频/GHz	1.2
波长/m	0.25
		阵元间距/m	0.125
阵元数	8
		脉冲重复频率/Hz	3000
发射带宽/MHz	2
		载机高度/km	5
脉冲数	8
		杂噪比	40
采样率/MHz	2
		载机速度/m/s	150
载机的偏航角/。	0

2.仿真结果：

图4表示天线阵元存在误差时，本发明与基于知识广义内积方法STAP输出信杂噪比结果图，图4横坐标表示归一化多普勒频率，纵坐标表示信杂噪比损失。图中的两条曲线分别表示本发明与广义内积方法STAP输出的信杂噪比示意图。

由图4中两种方法形成的信杂噪比损失可见，在天线阵元存在误差时的情况下，广义内积方法性能严重下降，这是因为根据理想杂波模型重构的杂波协方差矩阵与实际数据不匹配，导致其效果较差。而本发明通过对天线阵元误差进行估计，修正了按照理想模型得到的杂波协方差矩阵，有效提升了基于知识辅助广义内积方法的性能，改善了杂波抑制效果。

Claims (2)

1.一种基于杂波子空间的阵元误差估计方法，包括如下步骤：

(1)计算杂波空时导向矢量矩阵：

利用雷达先验信息，计算杂波空时导向矢量矩阵；

(2)奇异值分解：

2a)对雷达接收数据进行奇异值分解，得到接收数据左奇异值矢量矩阵，主、次奇异值对应的对角矩阵和接收数据右奇异值矢量矩阵；

2b)对杂波空时导向矢量矩阵进行奇异值分解，得到左奇异值矢量矩阵，主、次奇异值对应的对角矩阵和右奇异值矢量矩阵；

(3)计算主奇异值的数目：

按照下式，计算得到主奇异值的数目：

$r=\min (\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{\mathrm{\σ}}_i\≥0.999\underset{i=1}{\overset{\mathrm{NM}}{\mathrm{\Σ}}}{\mathrm{\σ}}_i)$

其中，r表示主奇异值的数目，m表示奇异值的数目，m的取值范围为：m＝1,…,NM，N表示天线阵元数，M表示雷达一个相干处理间隔内发射的脉冲数，i表示奇异值个数的标号，i的取值范围为：i＝1,…,m，σ_i表示对杂波空时导向矢量矩阵进行奇异值分解得到的奇异值，min(·)表示取最小值操作；

(4)构造正交补空间矩阵：

取出左奇异矢量矩阵，去掉左奇异矢量矩阵中主奇异值对应矩阵的列，将剩余列构成杂波的正交补空间矩阵；

(5)估计天线阵元误差：

5a)取出接收数据左奇异矢量矩阵中的首列，利用其与杂波的正交补空间矩阵的正交性，按照下式，估计天线阵元锥削矢量：

$\underset{\hat{t}}{\min }\left({\left|\right|U^H(e_1\·/\hat{t})\left|\right|}_2^2\right)$

其中，

表示估计的天线阵元锥削矢量，min(·)表示取最小值操作，U表示杂波的正交补空间矩阵，H表示矩阵共轭转置操作，e₁表示接收数据左奇异矢量矩阵中的首列，||·||₂表示2范数操作，·/表示对向量对应元素相除的点除操作；

5b)按照下式，估计天线阵元误差：

$\mathrm{\Δ}\hat{t}=\hat{t}(1:N)-1_N$

其中，表示估计的天线阵元误差，

表示估计的天线阵元锥削矢量，N表示天线阵元数，1_N表示N个全为1的列矢量；

(6)杂波抑制：

6a)利用M个全为1的列矢量构成的理想时域导向矢量，与估计出的天线阵元误差进行Kronecker积运算，得到修正后的空时锥削矢量，利用此空时锥削矢量与杂波空时导向矢量进行Hadamard积运算，得到修正后的杂波信号模型；

6b)利用杂波信号模型建立杂波协方差矩阵，进行杂波抑制。

2.根据权利要求1所述的基于杂波子空间的阵元误差估计方法，其特征在于：步骤（1）所述的计算杂波空时导向矢量矩阵，包括以下几个步骤：

第一步，按照下式，计算归一化多普勒频率：

其中，w₁表示归一化多普勒频率，v表示雷达载机速度，θ表示地面杂波散射体的俯仰角，

表示地面杂波散射体的方位角，λ表示雷达工作波长，f_r表示雷达脉冲发射的重复频率；

第二步，按照下式，计算归一化空间频率：

其中，w₂表示归一化空间频率，d表示天线阵元间距，θ表示地面杂波散射体的俯仰角，

表示地面杂波散射体的方位角，ψ表示载机速度方向与雷达天线阵列轴向的夹角，λ表示雷达工作波长；

第三步，按照下式，计算时域导向矢量：

$u\left(w_1\right)={\left[\begin{array}{ccccccc}1& e^{j2\mathrm{\π}{w}_1*1}& e^{j2\mathrm{\π}{w}_1*2}& \·\·\·& e^{j2\mathrm{\π}{w}_1*(m_1-1)}& \·\·\·& e^{j2\mathrm{\π}{w}_1(M-1)}\end{array}\right]}^T$

其中，u(w₁)表示时域导向矢量，w₁表示归一化多普勒频率，m₁表示归一化多普勒频率的数目，m₁的取值范围为：m₁＝1,…,M，M表示雷达一个相干处理间隔内发射的脉冲数，

表示第m₁个归一化多普勒频率对应的相位，T表示转置操作；

第四步，按照下式，计算空域导向矢量：

$u\left(w_2\right)={\left[\begin{array}{ccccccc}1& e^{j2\mathrm{\π}{w}_2*1}& e^{j2\mathrm{\π}{w}_2*2}& \·\·\·& e^{j2\mathrm{\π}{w}_2*(m_2-1)}& \·\·\·& e^{j2\mathrm{\π}{w}_2*(N-1)}\end{array}\right]}^T$

其中，u(w₂)表示空域导向矢量，w₂表示归一化空间频率，m₂表示归一化多普勒频率的数目，m₂的取值范围为：m₂＝1,…,N，

表示第m₂个归一化空间频率对应的相位，N表示雷达天线阵元数；

第五步，按照下式，计算空时导向矢量矩阵：

$V=u\left(w_1\right)\⊗u\left(w_2\right)$

其中，V表示空时导向矢量矩阵，u(w₁)表示时域导向矢量，w₁表示归一化多普勒频率，u(w₂)表示空域导向矢量，w₂表示归一化空间频率，

表示Kronecker积。

Images