一种冲击噪声下MIMO雷达阵列的诊断方法与流程

本发明涉及多输入多输出(mimo)雷达阵列诊断方法，具体涉及一种冲击噪声下mimo雷达阵列的诊断方法。

背景技术：

mimo(multiple-inputmultiple-output)雷达的诞生得益于mimo技术在无线通信中发展与延伸。相比于传统相控阵雷达，mimo雷达根据波形分集思想，由发射端发射正交波形，在接收端对接收信号进行匹配滤波处理，形成大孔径的虚拟阵列，以提高目标探测精度和参数估计能力。天线阵列作为雷达的重要组成部分，能实现自适应波束形成、空间谱估计和信源分离等功能，但受恶劣的天气、元器件腐蚀磨损、性能老化等因素的影响，使得雷达阵列中的收发组件出现故障，导致mimo雷达接收阵元无法准确接收目标信号，进而导致mimo雷达探测性能明显下降。由于对mimo雷达阵列中故障阵元进行维修的前提是需要确定其位置，因此，研究有效的雷达阵列诊断方法来确定故障阵元在阵列中的位置具有重要理论意义和实际应用价值。

针对于故障阵元位置诊断的问题，jiang等人在论文“impairedarraydiagnosisandmitigationwithkhatri–raoprocessing”(ieeeantennasandwirelesspropagationletters,2018,17(12):2354-2358)中提出一种对接收数据进行差分处理的均匀线性阵列诊断方法，利用故障阵元接收的数据与无故障阵元接收信号的相关性计算协方差矩阵，并向量化以构成虚拟阵列并经差分处理，利用处理后的数据诊断出故障阵元位置。zhuo等人在论文“arraydiagnosisandangleestimationinbistaticmimoradarunderarrayantennafailures”(ieeeradar,sonar&navigation,2019,13(7):1180-1188)中提出了一种基于图像熵的mimo雷达故障接收阵元诊断方法，该方法根据故障阵元和正常阵元对应的输出信号所包含信息量的不同，计算出各接收阵元对应的图像熵值，根据故障阵元与无故障阵元所对应图像块熵值的差异特性，以诊断出接收阵列中故障阵元的位置。zhang等人在论文“arraydiagnosisusingsignalsubspaceclusteringinmimoradar”(ieeeelectronicsletters,2020,56(2):99-102)中提出一种基于子空间聚类的mimo雷达故障阵元诊断方法，该方法先将虚拟阵列协方差矩阵经特征值分解以获得信号子空间，根据信号子空间各数据点的密度分布与距离特性采用密度峰值聚类算法，选取聚类中心，并划分出不同簇类，再确定异常簇类，最后，根据异常簇类确定故障阵元的位置。

以上方法对于mimo雷达的故障阵元诊断问题的研究大多数是以高斯噪声为假设前提。而雷达实际工作环境中，所遇到的噪声通常具有一定的冲击性，例如，环境噪声、大气噪声、雷达杂波、同信道射频干扰、低压电力线内部噪声以及系统电路上的瞬间尖峰等。冲击噪声(impulsivenoise)具有尖峰或厚重拖尾，通常用α稳定分布来表示，由于α稳定分布具有脉冲尖峰以及不具备二阶统计特性，导致现有基于二阶统计特性的jiang方法、zhang方法以及计算图像熵的zhuo方法性能下降甚至失效。

技术实现要素：

发明目的：本发明的目的在于提供一种能够准确诊断出mimo雷达故障阵元的位置的冲击噪声下mimo雷达阵列的诊断方法。

技术方案：本发明的冲击噪声下mimo雷达阵列的诊断方法，包括以下步骤：

(1)当mimo雷达阵列中存在故障阵元时，在mimo雷达接收阵列端获取接收数据矩阵；

(2)利用无穷范数归一化加权系数对存在故障阵元的mimo雷达接收信号进行预处理，得到加权处理后的mimo雷达接收阵列输出数据矩阵；

(3)对预处理后的mimo雷达接收信号进行匹配滤波，得到虚拟阵列在多个脉冲周期内的输出信号；

(4)对虚拟阵列输出信号矩阵进行取模处理，获得取模后的数据矩阵，并对取模后的数据矩阵量化处理得到数据矩阵，再将数据矩阵按序分割成多个图像块；

(5)对数据矩阵进行分析，依次计算各个接收与发射阵元对应的互信息值的均值，并构成接收阵列与发射阵列的互信息值矢量；

(6)选取合适的门限值对接收阵列互信息值矢量和发射阵列互信息值矢量分别进行门限检测，得到在矢量中低于门限值的元素及其在矢量中位置，根据得到的位置分别来确定接收和发射阵列中故障阵元的位置。

步骤(1)中，所述mimo雷达阵列接收的数据矩阵的计算公式为：

式中，表示n个接收阵元在第k个脉冲周期内的第p个采样数据；βl表示第l个目标的反射系数；b表示各发射阵元发射同频且相互正交的周期性编码信号矩阵；表示含故障阵元的接收阵列导向矢量；表示故障阵元的发射阵列导向矢量；表示冲击噪声矩阵。

步骤(2)中，所述无穷范数归一化加权系数的计算公式为：

式中，||·||∞表示无穷范数；表示矩阵的第n行第p列元素；|·|表示求模运算。

步骤(3)中，对预处理后的mimo雷达接收信号进行匹配滤波，获得m×n个虚拟阵元在k个脉冲周期的输出信号矩阵。

所述虚拟阵元输出信号矩阵的计算公式为：

式中，表示虚拟阵列输出信号矩阵；和分别表示含故障阵元的发射和接收阵列流型矩阵；s表示目标系数矩阵；表示噪声矩阵；表示khatri-rao积。

步骤(4)中，所述将数据矩阵按序分割成多个图像块具体为，将数据矩阵自上而下每m行分成一个子图像块，共形成n个子图像块矩阵；再将数据矩阵重排成新的数据矩阵，将自上而下每n行分成一个图像块，共构成m个图像块。

步骤(5)中，计算第n0个接收阵元图像块与其余接收阵元图像块的互信息值并取均值

所述接收阵元对应的互信息值的均值计算公式为：

式中，表示第n0(n0＝1,2,...,n)个接收阵元图像块与第n1(n1＝1,2,...,n,n1≠n0)个接收阵元图像块的互信息值。

步骤(5)中，计算第m0个发射阵元图像块与其余发射阵元图像块的互信息值并取均值

步骤(5)中，所述发射阵元对应的互信息值的均值计算公式为：

式中，表示第m0(m0＝1,2,...,m)个接收阵元图像块与第m1(m1＝1,2,...,n,m1≠m0)个接收阵元图像块的互信息值。

有益效果：本发明与现有技术相比，其有益效果在于：(1)利用了mimo雷达各阵元图像数据间相互包含的信息度量，充分挖掘各阵元图像数据之间的统计相关特性，从而诊断出mimo雷达阵列中正常阵元与故障阵元。相比于现有计算熵值的阵列诊断方法，本发明方法在低信噪比下的故障阵元诊断正确率更高，而且在较大范围的故障阵元数目下诊断性能具有较好的稳健性；(2)无需借助额外的测量设备来获取阵列远场或者近场方向图，仅需mimo雷达的接收数据即可完成发射和接收阵列中故障阵元位置的诊断，使得雷达具备故障阵元的自诊断能力，具有较好的实用性；(3)能够冲击噪声环境下准确诊断出mimo雷达故障阵元的位置，在低信噪比的冲击噪声环境下仍能保持较好的性能，从而有效解决了在mimo雷达在冲击噪声环境下故障阵元诊断问题。

附图说明

图1为本发明的流程示意图；

图2为本发明中双基地mimo雷达的结构图；

图3为本发明中无穷范数处理后的正常阵元和故障阵元的熵值和互信息值计算结果图；

图4为本发明中当mimo雷达阵列中含有2个位置随机出现的故障发射阵元和4个随机的故障接收阵元时故障阵元诊断成功概率随信噪比变化关系图；

图5为本发明中故障阵元诊断成功概率随故障阵元数的变化关系图。

具体实施方式

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

如图1所示，本发明的冲击噪声下mimo雷达阵列的诊断方法，当mimo雷达阵列中存在故障阵元时，故障阵元对应的图像数据中无目标信息，即故障阵元对应的图像数据中不包含与正常阵元图像数据中相似的目标信息，然而，正常阵元图像数据中包含与其余正常阵元图像数据相似的目标信息，因此，本发明采用各阵元对应的图像数据之间的互信息值来衡量各阵元所包含目标信息的差异程度，来进行故障阵元位置的诊断。具体包括以下步骤：

(1)当mimo雷达阵列中存在故障阵元时，mimo雷达接收阵列获取在第k个脉冲周期内的p个采样数据矩阵：

式中，表示n个接收阵元在第k个脉冲周期内的第p个采样数据；βl表示第l个目标的反射系数；b表示各发射阵元发射同频且相互正交的周期性编码信号矩阵；表示含故障阵元的接收阵列导向矢量；表示故障阵元的发射阵列导向矢量；表示冲击噪声矩阵；其中，qt为故障发射阵元位置集合，[at(φl)]m表示发射阵列导向矢量的第m(m＝1,2,…,m)个元素；[ar(θl)]n表示接收阵列导向矢量第n(n＝1,2,…,n)个元素，m和n分别为发射阵元和接收阵元数。

假设双基地mimo雷达系统有m个发射阵元和n个接收阵元，其收发阵列的阵元间距分别为dt和dr，且均为均匀线阵。在空间远场空间中有l个非相干观测目标，则第l(l＝1,...,l)个目标相对于发射阵列和接收阵列的波离方向角(dod)和波达方向角(doa)分别为(φl,θl)。各发射阵元发射同频且相互正交的周期性编码信号其中，bp＝[bp(1),…,bp(m)]^t，bp(m)为第m(m＝1,2,…,m)个发射阵元所发射的编码信号在第p(p＝1,2,…,p)个采样值，p为每个脉冲信号的采样个数，(·)^t表示转置。则mimo雷达接收阵列在第k个脉冲周期内的p个采样的输出信号为

式中，其中，dk(:,p)表示接收信号第k个脉冲周期内第p个采样数据；βl为第l个目标的反射系数；为接收阵列导向矢量；为发射阵列导向矢量，其中，λ表示发射信号波长；表示冲击噪声矩阵。

在实际工作环境中，mimo雷达天线阵列受恶劣环境以及阵列组件老化等因素的影响，会出现不同程度故障的阵元。当发射阵列中第qt(qt∈qt)个阵元出现故障时，则发射阵列导向矢量at(φl)中第qt个元素为零，即阵元故障发射阵列导向矢量的表达式为

式中，[at(φl)]m表示发射阵列导向矢量的第m个元素，qt表示故障发射阵元位置集合。同样，当接收阵列中第qr(qr∈qr)个阵元出现故障，接收阵列导向矢量ar(θl)中的第qr个元素为零，则阵元故障接收阵列导向矢量的表达式为，

式中，[ar(θl)]n表示接收阵列导向矢量第n个元素，qr表示故障接收阵元位置集合。则存在故障阵元时的mimo雷达在第k个脉冲周期内接收信号可表示为

(2)为了抑制冲击噪声的影响，利用无穷范数归一化加权系数对存在故障阵元的mimo雷达接收信号进行预处理，得到加权处理后的mimo雷达接收阵列输出信号矩阵，具体步骤如下：

对雷达接收数据用无穷范数进行归一化加权预处理，其归一化加权系数表达式为

式中，||·||∞表示无穷范数；表示矩阵的第n行第p列元素；|·|表示求模运算。第k个脉冲周期内p次采样的mimo雷达接收信号经过wk(p)加权处理后，可表示为

式中，为第l个目标的反射系数；为无穷范数归一化后的噪声数据。将式(6)所表示的p个采样值写成矩阵形式，即为

式中，为取第l个目标的反射系数构建对角矩阵，为零均值有限方差的噪声矩阵。

(3)对预处理后的mimo雷达接收信号进行匹配滤波获得m×n个虚拟阵元在k个脉冲周期的输出信号矩阵，具体步骤如下：

根据发射阵列各波形的正交性，用发射信号矩阵b对接收阵列信号进行匹配滤波可得虚拟阵列输出信号，即

式中，将按列写成mn×1维列矢量，即

式中，gk＝vec(gk)，vec(·)表示将矩阵向量化；表示khatri-rao积；和分别表示含故障阵元的发射和接收阵列流型矩阵；为一个n×l维矩阵,由l个接收阵列导向向量组成；为一个m×l维矩阵，由l个发射阵列导向向量组成；为第k个快拍的目标系数向量。虚拟阵列在k个快拍下的输出信号矩阵为

式中，为虚拟阵列输出信号矩阵；表示目标反射系数矩阵；g＝[g1,…,gk]，表示噪声矩阵。

(4)对虚拟阵列输出信号矩阵进行取模处理，获得取模后的数据矩阵，并对取模后的数据矩阵量化处理得到数据矩阵，再将数据矩阵按序分割成多个图像块；具体计算过程如下：

对虚拟阵列输出信号矩阵进行取模处理，获得求模后的数据矩阵即等效于灰度值为非整数连续量的灰度图像数据矩阵，对其进行量化处理，将其灰度值映射到0～255之间，得到离散的灰度图像数据：

式中，y表示离散的灰度图像数据矩阵；round(·)表示四舍五入运算。

当阵列中存在故障阵元时，一组故障阵元图像数据中包含的关于另一组阵元图像数据的信息度量较小，故可利用互信息值来衡量图像块之间相互包含的信息度量差异。为了进行接收阵列中故障阵元位置的诊断，将处理后的数据矩阵y自上而下每m行分成一个子图像块，共构成n个接收阵元图像块yn(n＝1,2,…,n)，即同样,为了进行发射阵列中故障阵元位置的诊断，将数据矩阵y重排成新的数据矩阵yt，即

yt＝[y1,ym+1,…,y(n-1)m+1,y2,ym+2,…,y(n-1)m+2,…,ym,y2m,…,ynm]^t(12)

式中，yi为矩阵y中第i(i＝1,2,…,mn)个行矢量。将yt并自上而下每n行分成一个图像块，共构成m个图像块，即

(5)对数据矩阵进行分析，依次计算各个接收与发射阵元对应的互信息值的均值，并构成接收阵列与发射阵列的互信息矢量；具体计算过程如下：

依次计算接收阵元图像块yn(n＝1,2,…,n)的熵值，即第n(n＝1,2,...,n)个接收阵元图像块yn的熵值为

式中，表示灰度值a∈[0,255]在图像块yn中出现的概率；fa表示灰度值a出现的次数。第n0(n0＝1,2,...,n)个接收阵元图像块与第n1(n1＝1,2,...,n,n1≠n0)个接收阵元图像块的互信息为

式中，表示第n0(n0＝1,2,...,n)个接收阵元图像块的熵值，表示第n1(n1＝1,2,...,n)个接收阵元图像块的熵值；则表示第n0(n0＝1,2,...,n)个接收阵元图像块与第n1(n1＝1,2,...,n,n1≠n0)个接收阵元图像块的联合熵，其表达式为

式中,pa,b为图像块中灰度值为a∈[0,255]和图像块中同名灰度值为b∈[0,255]时的联合概率，且满足为了准确诊断出故障接收阵元的位置，将第n0(n0＝1,...,n)个接收阵元图像块与其余接收阵元图像块的互信息值取均值，可表示为

同理，对数据矩阵yt进行分析，计算第m0(m0＝1,2,...,m)个发射阵元图像块与其余发射阵元图像块的互信息值并取均值，即

同样，将第m0(m0＝1,2,...,m)个发射阵元图像块与其余发射阵元图像块的互信息值取均值，即为

式中，表示第m0(m0＝1,2,...,m)个接收阵元图像块与第m1(m1＝1,2,...,n,m1≠m0)个接收阵元图像块的互信息值。

当接收阵列中第n0(n0∈qr)个接收阵元为故障阵元时，该阵元对应图像块与其余阵元图像块的互信息均值较小，当接收阵列中第个接收阵元为正常阵元时，其图像块与其余阵元图像块的互信息均值较大。同理当发射阵列中第m0(m0∈qt)个发射阵元位出现故障时，其发射阵元图像块对应的互信息均值较小，而当发射阵列中第个发射阵元为正常阵元时，其图像块对应的互信息均值较大。

接着，将得到的接收阵列与发射阵列中各阵元对应的互信息均值和构成互信息值矢量v＝[v1,v2,…,vn]和u＝[u1,u2,…,um]。

(6)选取合适的门限值对接收阵列互信息值矢量v和发射阵列互信息值矢量u分别进行门限检测，得到在矢量v和u中低于门限值的元素及其在矢量中位置，根据该位置来确定接收和发射阵列中故障阵元的位置。

将接收阵列与发射阵列中各阵元图像块的互信息均值和分别组成接收阵列的互信息矢量v＝[v1,v2,…,vn]和发射阵列的互信息值矢量u＝[u1,u2,…,um]。由于故障接收阵元对应的互信息均值要小于正常接收阵元对应的互信息均值，因此，设置合适的门限值vr＝β·mean(v)对矢量v进行门限检测，其中，β表示调节因子，mean(·)表示取均值。当矢量v中的值小于阈值vr时，则该值对应的下标为故障接收阵元在接收阵列中的位置。同理可知，矢量u中小于门限值ut＝β·mean(u)的值对应的下标即为故障发射阵元在发射阵列中的位置。最后可获得故障接收阵元和故障发射阵元的位置集合估计值和最终实现mimo雷达发射和接收阵列中故障阵元位置的诊断。

本发明的技术效果可以通过以下仿真结果进一步说明。为了验证本发明在冲击噪声下的故障阵元诊断的性能，设计了熵值方法、子空间聚类算法、差分处理方法以及无穷范数熵值方法与无穷范数互信息方法的对比实验。如图2所示，假设双基地mimo雷达发射阵元数m＝10，接收阵元数n＝20，目标个数l＝3，各发射阵列的发射角φ＝{20°，-25°,30°},各接收阵列的接收角θ＝{-13°，-5°,-15°}，每个周期内采样个数为p＝256，无穷范数互信息方法中用于门限检测的调节因子β＝0.98。仿真实验中的噪声属于冲击噪声，冲击噪声的分布特性是由特征指数α和分散系数γ确定，当α＜2时，其方差是无穷大的，因此，可以用广义信噪比(gsnr)来代替一般信噪比(snr)，即信号平均功率与分散系数γ的比值，其表达式如下：

式中，||·||f表示frobenius范数。定义故障阵元的诊断准确率为

式中，ct为故障发射阵元诊断成功的次数；ct为故障接收阵元诊断成功的次数；t为蒙特卡洛实验次数。

仿真实验1：如图3所示，假设接收阵列中位置{3,6,9,12,17}的阵元为故障阵元，广义信噪比gsnr＝10db，脉冲周期数k＝200，特征指数α＝1.6。当故障阵元对应的图像数据中包含冲击噪声的脉冲尖峰时，其信息量的不确定性增大，导致该故障阵元所包含的信息度量变大，致使熵值方法和互信息方法无法利用信息度量来正确区分正常阵元和故障阵元，因此，需要对mimo雷达接收数据进行无穷范数处理以抑制冲击噪声。由图3可知，位置{3,6,9,12,17}的故障接收阵元图像数据的熵值和互信息值分别低于其他正常接收阵元图像数据的熵值和互信息值，值得注意的是，采用互信息对各阵元图像数据进行信息度量时，故障阵元与正常阵元的特征区分度明显大于熵值度量方法。这是由于熵值方法仅仅计算了各阵元图像数据的信息度量，而互信息方法计算了各阵元图像数据间相互包含的信息度量，即分析了各阵元图像数据之间的统计相关程度。因此，相比于熵值方法，互信息方法更利于进行故障阵元位置的诊断。

仿真实验2：如图4所示，假设mimo雷达阵列中含有2个位置随机出现的故障发射阵元和4个随机的故障接收阵元，广义信噪比gsnr的变化范围为0～20db，蒙特卡洛实验次数t＝100，脉冲周期数k＝200。图4(a)中，冲击噪声特征指数α＝1.6，图4(b)中，冲击噪声特征指数α＝1.2。由图4可知，由于冲击噪声不具备二阶统计量的特性，导致差分处理方法和子空间聚类算法在不同强度的冲击噪声下对故障阵元诊断均失效。对于熵值方法，当冲击噪声特征指数α＝1.6时，在信噪比低于11db时，该方法的诊断准确率会逐渐下降，当冲击噪声特征指数α＝1.6变为α＝1.2时，即冲击噪声的脉冲尖峰幅度增大，使得故障阵元所包含的信息度量变大，从而导致熵值方法诊断失效。在不同的α值下，随着信噪比的不断增加，无穷范数互信息方法和无穷范数熵值方法的故障阵元诊断正确率均能不断提高，尤其在低信噪比的冲击噪声环境下，无穷范数互信息方法的诊断性能要优于无穷范数熵值方法。

仿真实验3：为了进一步验证无穷范数互信息方法在不同故障阵元数时的诊断稳健性，假设接收阵元中位置随机的故障阵元个数由1～12依次增加，广义信噪比gsnr＝7db，冲击噪声特征指数α＝1.6，其余参数与实验2一致。如图5所示，熵值方法、子空间聚类算法及差分处理方法在不同故障阵元数下完全失效。无穷范数熵值方法当故障阵元数超过8个时，诊断性能急剧恶化直至失效，但是无穷范数互信息方法在故障阵元数达到11个时仍能保持100％的诊断准确率，因此无穷范数互信息方法具有良好的稳健性。