CN105869189A - 基于ffdiag算法的雷达目标盲源分离方法 - Google Patents

Abstract

本发明公开了一种基于FFDIAG算法的雷达目标盲源分离方法，其主要思路为：依次计算基于F‑范数的目标矩阵代价函数L(V)和第k个雷达目标矩阵经过第n次迭代后得到的N×N维解混迭矩阵V_(n+1)，并将V_(n+1)写成与W_(n)有关的代价函数形式L(W_(n))，进而计算第n次迭代后包含的代价函数的绝对值；若所述绝对值小于等于ε，依次计算第n次迭代后第k个雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵W_(n)、第k个雷达目标矩阵经过第n次迭代后得到的N×N维解混迭矩阵V_(n+1)和第k个雷达目标矩阵经过第n次迭代后的雷达目标矩阵对角化结构并将作为第k个雷达目标矩阵的联合对角化结构，进而得到K个雷达目标矩阵各自对应的联合对角化结构并解决K个雷达目标矩阵中的盲源分离。

Description

基于FFDIAG算法的雷达目标盲源分离方法

技术领域

本发明属于基于联合对角化解盲源分离领域，涉及一种基于FFDIAG算法的雷达目标盲源分离方法，即基于快速Frobenius联合近似对角化算法的快速复数域联合对角化方法，适用于求解复数域的盲源分离问题。

背景技术

现有的联合对角化算法大多都是规定雷达目标矩阵组中必需含有一个正定矩阵，这样就可以使用正定矩阵对雷达目标信号进行预白化处理，从而将待求的联合对角化器转换成正交(酉)矩阵并进行求解，如JADE算法和JMMD算法；但是由于实际情况中的雷达目标矩阵通常都是由统计方式获得，因而具有不同程度的误差，并且正定性不能得以保证，进而造成白化处理的不准确，另外在预白化处理过程中带入的误差也很难在后续算法中加以校正，使得总体算法性能随之大打折扣。于是，许多学者提出了非正交联合对角化算法，该非正交联合对角化算法避免了白化预操作，如基于Log-Likelihood的近似联合对角化算法和SVDJD算法等，但该非正交联合对角化算法对雷达目标矩阵仍有约束，且规定雷达目标矩阵均为正定俄米特矩阵，而在盲源分离的情况中，该约束也很难满足，从而限定了非正交联合对角化算法的适用领域。

近几年来，又有研究者提出了非正交联合对角化的改进方法，不仅避免了预白化处理，而且不再要求雷达目标矩阵的正定性，从而拓展了非正交联合对角化算法的应用范围，其中相对有名的算法包括子空间拟合算法、ACDC算法、FFDIAG算法、QRJ2D算法、LSA(B)算法、J-Di算法以及DNJD算法；然而，在此类算法中，大多数都是假定雷达目标矩阵组或者(解)混迭矩阵为实值，而在实际应用中数据形式往往是复值数据的，如雷达信号处理或生物医学信号；另外，在一些具体问题如频率估计或盲波束形成中，不可避免地需要对复值数据进行处理。

此外，在众多联合对角化方法中，并不是每一个方法都能够从实数域延伸到复数域，因为联合对角化方法的核心步骤是假设基于实数，若此假设不存在，会使得对应算法的后续步骤无法展开；例如，Ziehe等研究人员提出的FFDIAG算法为基于实数假设的算法，该算法不仅避免了预白化处理步骤，也不要求雷达目标矩阵组以及(解)混迭矩阵的正定性，是一个简单高效的算法；但是该算法并不能够在复数域的盲源分离情况中直接使用，并且关键步骤的成立是以实数假设为前提，因此极大的限制了该算法的应用价值。

发明内容

针对以上现有技术存在的不足，本发明的目的在于提出一种基于FFDIAG算法的雷达目标盲源分离方法，该种基于FFDIAG算法的雷达目标盲源分离方法是一种适用于复数域的联合对角化方法，计算复杂度低，不仅能够避免预白化处理和正定性限制，而且将适用范围从实数域扩展到复数域，是对FFDIAG算法的一种传承，因而命名为CVFFDIAG(Complex-Valued FFDIAG)算法。

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

一种基于FFDIAG算法的雷达目标盲源分离方法，包括以下步骤：

步骤1，分别设定雷达目标矩阵总个数K以及N×N维复矩阵集合C^N×N，并计算得到基于F-范数的雷达目标矩阵代价函数L(V)；其中，V表示设定的解混迭矩阵，K和N分别为自然数；

步骤2，初始化：令n为迭代次数，且初始值为1；k表示第k个雷达目标矩阵，且k的初始值为1，V₍₁₎表示N×N维单位阵，W₍₀₎表示N×N维零矩阵；

步骤3，计算得到第k个雷达目标矩阵经过第n次迭代后得到的N×N维解混迭矩阵V_(n+1)，并将V_(n+1)写成与W_(n)有关的代价函数形式L(W_(n))；其中，W_(n)表示第n次迭代后第k个雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵，N为自然数；

步骤4，根据与W_(n)有关的代价函数形式L(W_(n))，计算得到关于的代价函数其中，令 表示w_ij的实部，表示w_ji的实部，表示w_ij的虚部，表示w_ji的虚部，w_ij表示W_(n)的第i行、第j列元素，w_ji表示W_(n)的第j行、第i列元素，W_(n)表示第n次迭代后第k个雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵，上标T表示转置；

步骤5，将关于的代价函数对求导，并令导数等于零，得到的表达式然后令i的初值为1，j∈{i+1,…,N}，得到的N-1个对应值；再令i加1，j∈{i+1,…,N}，得到的N-2个对应值；直到i＝N-1，j∈{i+1,…,N}，得到的1个对应值；

利用i∈{1,…,N-1}、j∈{i+1,…,N}时对应的所有值，计算得到第n次迭代后包含的代价函数的绝对值并设定收敛条件：

$\left|\min L\left({\hat{w}}_{ij}\right)\right|\≤\mathrm{\ϵ}$

其中，ε表示设定的趋于无穷小的极小数；将所述第n次迭代后包含的代价函数的绝对值与ε进行比较；若所述第n次迭代后包含的代价函数的绝对值大于ε，则令n加1，返回步骤3；

若所述第n次迭代后包含的代价函数的绝对值小于等于ε，则依次计算得到第n次迭代后第k个雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵W_(n)、第k个雷达目标矩阵经过第n次迭代后得到的N×N维解混迭矩阵V_(n+1)和第k个雷达目标矩阵经过第n次迭代后的雷达目标矩阵对角化结构并将所述第k个雷达目标矩阵经过第n次迭代后的雷达目标矩阵对角化结构作为第k个雷达目标矩阵的联合对角化结构；然后令k加1，同时初始化n为1，返回步骤3；直到k＝K，进而得到K个雷达目标矩阵各自对应的联合对角化结构，用于解决K个雷达目标矩阵中的盲源分离问题。

本发明方法相对于现有算法的优点在于：

第一，本发明的创新点在于将一种联合对角化算法从实数域延伸到了复数域；

第二，本发明方法能够避免对雷达雷达目标矩阵及(解)混迭矩阵的约束，产生了更广泛的适用范围和更大的实用价值。

第三，本发明方法免除了预白化处理和正定性约束，减小了运算量，节省了时间。

附图说明

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

图1是本发明的一种基于FFDIAG算法的雷达目标盲源分离方法流程图；

图2A是未加入噪声的全局拒绝水平随迭代次数变化的曲线示意图，

图2B是未加入噪声的对角化误差随迭代次数变化的曲线示意图；

图3A是加入噪声的全局拒绝水平随迭代次数变化的曲线示意图，

图3B是加入噪声的平均全局拒绝水平随迭代次数变化的曲线示意图；

图4是分别使用DC算法、FAJD算法和本发明方法得到的全局拒绝水平随迭代次数变化的曲线示意图。

具体实施方式

参照图1，为本发明的一种基于FFDIAG算法的雷达目标盲源分离方法示意图，该种基于FFDIAG算法的雷达目标盲源分离方法，包括以下步骤：

步骤1，分别设定雷达雷达目标矩阵总个数K以及N×N维复矩阵集合C^N×N，并计算得到基于F-范数的雷达目标矩阵代价函数L(V)；其中，V表示设定的解混迭矩阵，K和N分别为自然数。

具体地，分别设定雷达目标矩阵总个数K以及N×N维复矩阵集合C^N×N，并计算得到基于F-范数的雷达目标矩阵代价函数L(V)，其表达式为：

$L\left(V\right)=\underset{V\∈C^{N\×N}}{min}\underset{k=1}{\overset{K}{\Σ}}off\left({\mathrm{VC}}^k{V}^H\right)---\left(1\right)$

其中，min表示求取最小值操作，C^N×N表示N×N维复矩阵集合，k∈{1,…,K}，K表示雷达目标矩阵总个数，V表示设定的解混迭矩阵，off(·)表示矩阵中所有非对角线元素的F-范数之和，C^k表示第k个雷达目标矩阵，上标H表示共轭转置，K和N分别为自然数。

步骤2，初始化：令n为迭代次数，且初始值为1；k表示第k个雷达目标矩阵，且k的初始值为1，V₍₁₎表示N×N维单位阵，W₍₀₎表示N×N维零矩阵。

步骤3，计算得到第k个雷达目标矩阵经过第n次迭代后得到的N×N维解混迭矩阵V_(n+1)，并将V_(n+1)写成与W_(n)有关的代价函数形式L(W_(n))；其中，W_(n)表示第n次迭代后第k个雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵，N为自然数。

具体地，本发明方法通过衡量雷达目标矩阵与对角矩阵的偏离程度求解最优解，即采用非对角线元素的平方和表征联合对角化的近似程度，所述基于F-范数的雷达目标矩阵代价函数L(V)用以求解设定的解混迭矩阵V，并通过乘性迭代机制对设定的解混迭矩阵V进行更新，计算得到第k个雷达目标矩阵经过第n次迭代后得到的N×N维解混迭矩阵V_(n+1)，其表达式为：

V_(n+1)＝(I+W_(n))V_(n)(2)

其中，I表示N×N维单位阵，下标n为迭代次数，且n的初值为1，V_(n)表示第k个雷达目标矩阵经过第n-1次迭代后得到的N×N维解混迭矩阵，V₍₁₎表示N×N维单位阵，W₍₀₎表示N×N维零矩阵；W_(n)表示第n次迭代后第k个雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵，V_(n+1)表示第k个雷达目标矩阵经过第n次迭代后得到的N×N维解混迭矩阵。

定义k∈{1,2,…,K}，表示第k个雷达目标矩阵经过第n次迭代后的雷达目标矩阵，C^k表示第k个雷达目标矩阵，V_(n+1)表示第k个雷达目标矩阵经过第n次迭代后得到的N×N维解混迭矩阵，上标H表示共轭转置；其中，第k个雷达目标矩阵的初值

对进行进一步代换近似，进而将所述基于F-范数的雷达目标矩阵代价函数L(V)变换为与W_(n)有关的代价函数形式L(W_(n))，W_(n)表示第n次迭代后第k个雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵。

在迭代过程中，通过对V_(n)可逆性的保证避免式(1)收敛到平凡解，即V＝0，必须保证V_(n)的可逆性，V_(n)表示第k个雷达目标矩阵经过第n-1次迭代后得到的N×N维解混迭矩阵，下面根据定义1和定理1确保I+W_(n)的可逆性。

定义1：一个N×N维矩阵F，若F所有对角线元素满足条件则F是严格对角占优的。

定理1(Levi-Desplanques理论)：如果一个N×N维矩阵F是严格对角占优的，则F是可逆矩阵。

根据定理1，如果I+W_(n)是对角占优的，则I+W_(n)是可逆的，又因为I+W_(n)所有对角线元素为1，那么只要满足则I+W_(n)满足对角占优，其中表示第n次迭代后雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵W_(n)中第i行、第j列元素；如果设定一个常数θ，0<θ<1，并令或者直接令使得I+W_(n)对角占优；此处设定θ＝0.9，||·||_F表示F范数，表示逼近。

其中，第k个雷达目标矩阵的初值在不失一般性的前提下，以第n步迭代为例进行阐述；在第n步迭代中，通过使用W_(n)来表示V_(n+1)，从而将第k个雷达目标矩阵经过第n-1次迭代后的雷达目标矩阵变换为更具对角化结构的表示第k个雷达目标矩阵经过第n次迭代后的雷达目标矩阵，表示用W_(n)替换V_(n+1)的第k个雷达目标矩阵经过n次迭代后的雷达目标矩阵，其表达式为：

$\stackrel{\&OverBar;}{C_{(n+1)}^k}=(I+W_{\left(n\right)})C_{\left(n\right)}^k{(I+W_{\left(n\right)})}^H---\left(3\right)$

令的对角部分为非对角部分为得到用和替换的第k个雷达目标矩阵经过n次迭代后的雷达目标矩阵其表达式为：

$\stackrel{\&OverBar;}{\stackrel{\&OverBar;}{C_{(n+1)}^k}}=D_{\left(n\right)}^k+D_{\left(n\right)}^k{W}_{\left(n\right)}^H+E_{\left(n\right)}^k+E_{\left(n\right)}^k{W}_{\left(n\right)}^H+W_{\left(n\right)}{D}_{\left(n\right)}^k+W_{\left(n\right)}{D}_{\left(n\right)}^k{W}_{\left(n\right)}^H+W_{\left(n\right)}{E}_{\left(n\right)}^k+W_{\left(n\right)}{E}_{\left(n\right)}^k{W}_{\left(n\right)}^H---\left(4\right)$

其中，表示第k个雷达目标矩阵经过第n-1次迭代后的雷达目标矩阵；此处，假设W_(n)和分别具有很小的模值，因而忽略上式中含二次W_(n)和的项，得到分别忽略含W_(n)和二次项的第k个雷达目标矩阵经过n次迭代后的雷达目标矩阵其表达式为：

由于式(5)中的为对角矩阵，因此只要项更加接近对角化结构，也就越接近对角化结构；即将第n次迭代后的N×N维解混迭矩阵V_(n+1)写成与W_(n)有关的代价函数形式L(W_(n))，其表达式为：

其中，I表示N×N维单位阵，V_(n)表示第k个雷达目标矩阵经过第n-1次迭代后得到的N×N维解混迭矩阵，off(·)表示矩阵中所有非对角线元素的F-范数之和，min表示求取最小值操作，C^N×N表示N×N维复矩阵集合，W_(n)表示第n次迭代后第k个雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵，表示的对角部分，表示的非对角部分，表示第k个雷达目标矩阵经过第n-1次迭代后的雷达目标矩阵，上标H表示共轭转置。

步骤4，根据与W_(n)有关的代价函数形式L(W_(n))，计算得到关于的代价函数其中，令 表示w_ij的实部，表示w_ji的实部，表示w_ij的虚部，表示w_ji的虚部，w_ij表示W_(n)的第i行、第j列元素，w_ji表示W_(n)的第j行、第i列元素，W_(n)表示第n次迭代后第k个雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵，上标T表示转置。

具体地，根据与W_(n)有关的代价函数形式L(W_(n))，计算得到关于的代价函数其过程为：

首先，将式(6)用元素的形式表达出来，得到L(W_(n))的元素表现形式

$L\left(\widetilde{W_n}\right)=\underset{k-1}{\overset{K}{\&Sigma;}}\underset{i=1}{\overset{N-1}{\&Sigma;}}\underset{j=i+1}{\overset{N}{\&Sigma;}}\left\{\right|w_{ij}{d}_j^k+w_{ji}^{*}{d}_i^k+e_{ij}^k{|}^2+|w_{ji}{d}_i^k+w_{ij}^{*}{d}_j^k+e_{ji}^k{|}^2\}---\left(7\right)$

其中，w_ij表示W_(n)的第i行、第j列元素，表示E^k的第i行、第j列元素，表示D^k对角线上的第i个元素，w_ji表示W_(n)的第j行、第i列元素，表示E^k的第j行、第i列元素，表示D^k对角线上的第j个元素，W_(n)表示第k个雷达目标矩阵中对角线元素分别为零的N×N维矩阵，D^k表示第k个雷达目标矩阵C^k的对角部分，E^k表示第k个雷达目标矩阵C^k的非对角部分，W_(n)表示第n次迭代后第k个雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵。

然后把式(7)分解成只包含w_ij、w_ji的代价函数L(w_ij,w_ji)，并将w_ij、w_ji各自的实虚部分离，进而计算得到关于的代价函数其中，表示w_ij和w_ji各自对应的实部和虚部组成的向量，w_ij表示W_(n)的第i行、第j列元素，w_ji表示W_(n)的第j行、第i列元素，W_(n)表示第n次迭代后第k个雷达目标矩阵中对角线元素分别为零的N×N维矩阵。

进而将式(7)分解成N(N-1)/2个子代价函数，得到L(W_n)的子代价函数形式L(w_ij,w_ji)：

$L(w_{ij},w_{ji})=\underset{k=1}{\overset{K}{\&Sigma;}}|w_{ij}{d}_j^k+w_{ji}^{*}{d}_i^k+e_{ij}^k{|}^2+|w_{ji}{d}_i^k+w_{ij}^{*}{d}_j^k+e_{ji}^k{|}^2,$

i∈{1,…,N-1}，j∈{i+1,…,N} (8)

显然，式(8)中只包含w_ij、w_ji两个复数未知数，w_ij表示W_(n)的第i行、第j列元素，w_ji表示W_(n)的第j行、第i列元素，W_(n)表示第n次迭代后第k个雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵。

分别将w_ij、和各自对应的实虚部分离，即将w_ij、和各自对应的实部分别用 和表示，将w_ij、和各自对应的虚部分别用 和表示，得到L(w_ij,w_ji)的另一种表现形式

其中，i∈{1,…,N-1}，j∈{i+1,…,N}。

对式(9)使用向量形式进行简化，首先分别令：

其中，上标T表示转置，表示w_ij的实部，表示w_ji的实部，表示w_ij的虚部，表示w_ji的虚部，w_ij表示W_(n)的第i行、第j列元素，w_ji表示W_(n)的第j行、第i列元素，W_(n)表示第n次迭代后第k个雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵；表示的实部，表示的实部，表示的虚部，表示的虚部，表示D^k对角线上的第i个元素，表示D^k对角线上的第j个元素，D^k表示第k个雷达目标矩阵C^k的对角部分。

进而将式(9)转化为关于的代价函数

其中，令 表示w_ij的实部，表示w_ji的实部，表示w_ij的虚部，表示w_ji的虚部，w_ij表示W_(n)的第i行、第j列元素，w_ji表示W_(n)的第j行、第i列元素，W_(n)表示第n次迭代后第k个雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵，上标T表示转置。

步骤5，将关于的代价函数对求导，并令导数等于零，得到的表达式令i的初值为1，j∈{i+1,…,N}，得到的N-1个对应值；再令i加1，j∈{i+1,…,N}，得到的N-2个对应值；直到i＝N-1，j∈{i+1,…,N}，得到的的1个对应值。

利用i∈{1,...,N-1}、j∈{i+1,...,N}时对应的所有值，计算得到第n次迭代后包含的代价函数的绝对值并设定收敛条件：

$\left|\min L\left(\widehat{w_{ij}}\right)\right|\&le;\mathrm{\&epsiv;}$

其中，ε表示设定的趋于无穷小的极小数，本发明设定ε等于10^-6；将所述第n次迭代后包含的代价函数的绝对值与极小数ε进行比较；若所述第n次迭代后包含的代价函数的绝对值大于ε，算法未收敛，则令n加1，返回步骤3。

若所述第n次迭代后包含的代价函数的绝对值小于等于ε，则算法收敛，此时依次计算得到第n次迭代后第k个雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵W_(n)、第k个雷达目标矩阵经过第n次迭代后得到的N×N维解混迭矩阵V_(n+1)和第k个雷达目标矩阵经过第n次迭代后的雷达目标矩阵对角化结构并将所述第k个雷达目标矩阵经过第n次迭代后的雷达目标矩阵对角化结构作为第k个雷达目标矩阵的联合对角化结构；然后令k加1，同时初始化n为1，返回步骤3；直到k＝K，进而得到K个雷达目标矩阵各自对应的联合对角化结构，用于解决K个雷达目标矩阵中的盲源分离问题。

具体地，将关于的代价函数对进行求导，并令其导数等于零，得到：

其中， 上标T表示转置，表示w_ij的实部，表示w_ji的实部，表示w_ij的虚部，表示w_ji的虚部，w_ij表示W_(n)的第i行、第j列元素，w_ji表示W_(n)的第j行、第i列元素；表示的实部，W_(n)表示第n次迭代后第k个雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵；表示的实部，表示的虚部，表示的虚部，表示D^k对角线上的第i个元素，表示D^k对角线上的第j个元素，D^k表示第k个雷达目标矩阵C^k的对角部分。

为了简化形式令

$B=\underset{k=1}{\overset{K}{\&Sigma;}}\underset{p=1}{\overset{4}{\&Sigma;}}{d}_p\left(k\right)d_p^T\left(k\right)---\left(13\right)$

其中，K表示设定的雷达目标矩阵总个数，上标T表示转置；进而将式(10)转化为：

因此得到的表达式

$\stackrel{\mathrm{\&Delta;}}{w_{ij}}=B^{-1}c---\left(15\right)$

在计算的表达式时利用d_p(k),p∈{1,…,4}，k∈{1,…,K}的内部特点以及B的块对角结构，即为了易于表示，令那么：

${\tilde{d}}_1\left(k\right)={\&lsqb;x_1\left(k\right),x_2\left(k\right),-x_3\left(k\right),x_4\left(k\right)\&rsqb;}^T$

${\tilde{d}}_2\left(k\right)={\&lsqb;x_3\left(k\right),x_4\left(k\right),x_1\left(k\right),-x_2\left(k\right)\&rsqb;}^T$

${\tilde{d}}_3\left(k\right)={\&lsqb;x_1\left(k\right),x_2\left(k\right),x_3\left(k\right),-x_4\left(k\right)\&rsqb;}^T$

${\tilde{d}}_4\left(k\right)={\&lsqb;x_3\left(k\right),x_4\left(k\right),-x_1\left(k\right),x_2\left(k\right)\&rsqb;}^T---\left(16\right)$

再令

${\mathrm{\&alpha;}}_{qp}={\mathrm{\&alpha;}}_{pq}=\underset{k=1}{\overset{K}{\&Sigma;}}{x}_p\left(k\right)x_q\left(k\right)---\left(17\right)$

其中，p∈{1,…,4}，x_p(k)对应表示为x₁(k)、x₂(k)、x₃(k)或x₄(k)，q∈{1,…,4}；x_q(k)对应表示为x₁(k)、x₂(k)、x₃(k)或x₄(k)，α_pq表示x_p(k)和x_q(k)相乘，α_qp表示x_q(k)和x_p(k)相乘；进而得到B的矩阵表达式

$\begin{array}{c}\hat{B}=\underset{k=1}{\overset{K}{\&Sigma;}}\&lsqb;{\tilde{d}}_1\left(k\right){{\tilde{d}}_1}^T\left(k\right)+{\tilde{d}}_2\left(k\right){{\tilde{d}}_2}^T\left(k\right)+{\tilde{d}}_3\left(k\right){{\tilde{d}}_3}^T\left(k\right)+{\tilde{d}}_4\left(k\right){{\tilde{d}}_4}^T\left(k\right)\&rsqb;\\ =\left[\begin{array}{cccc}{\mathrm{\&alpha;}}_{11}& {\mathrm{\&alpha;}}_{12}& -{\mathrm{\&alpha;}}_{13}& {\mathrm{\&alpha;}}_{14}\\ {\mathrm{\&alpha;}}_{12}& {\mathrm{\&alpha;}}_{22}& -{\mathrm{\&alpha;}}_{23}& {\mathrm{\&alpha;}}_{24}\\ -{\mathrm{\&alpha;}}_{13}& -{\mathrm{\&alpha;}}_{23}& {\mathrm{\&alpha;}}_{33}& -{\mathrm{\&alpha;}}_{34}\\ {\mathrm{\&alpha;}}_{14}& {\mathrm{\&alpha;}}_{24}& -{\mathrm{\&alpha;}}_{34}& {\mathrm{\&alpha;}}_{44}\end{array}\right]+\left[\begin{array}{cccc}{\mathrm{\&alpha;}}_{33}& {\mathrm{\&alpha;}}_{34}& {\mathrm{\&alpha;}}_{13}& -{\mathrm{\&alpha;}}_{23}\\ {\mathrm{\&alpha;}}_{34}& {\mathrm{\&alpha;}}_{44}& {\mathrm{\&alpha;}}_{14}& -{\mathrm{\&alpha;}}_{24}\\ {\mathrm{\&alpha;}}_{13}& {\mathrm{\&alpha;}}_{14}& {\mathrm{\&alpha;}}_{11}& -{\mathrm{\&alpha;}}_{12}\\ -{\mathrm{\&alpha;}}_{23}& -{\mathrm{\&alpha;}}_{24}& -{\mathrm{\&alpha;}}_{12}& {\mathrm{\&alpha;}}_{22}\end{array}\right]\\ +\left[\begin{array}{cccc}{\mathrm{\&alpha;}}_{11}& {\mathrm{\&alpha;}}_{12}& {\mathrm{\&alpha;}}_{13}& -{\mathrm{\&alpha;}}_{14}\\ {\mathrm{\&alpha;}}_{12}& {\mathrm{\&alpha;}}_{22}& {\mathrm{\&alpha;}}_{23}& -{\mathrm{\&alpha;}}_{24}\\ {\mathrm{\&alpha;}}_{13}& {\mathrm{\&alpha;}}_{23}& {\mathrm{\&alpha;}}_{33}& -{\mathrm{\&alpha;}}_{34}\\ -{\mathrm{\&alpha;}}_{14}& -{\mathrm{\&alpha;}}_{24}& -{\mathrm{\&alpha;}}_{34}& {\mathrm{\&alpha;}}_{44}\end{array}\right]+\left[\begin{array}{cccc}{\mathrm{\&alpha;}}_{33}& {\mathrm{\&alpha;}}_{34}& -{\mathrm{\&alpha;}}_{13}& {\mathrm{\&alpha;}}_{23}\\ {\mathrm{\&alpha;}}_{34}& {\mathrm{\&alpha;}}_{44}& -{\mathrm{\&alpha;}}_{14}& {\mathrm{\&alpha;}}_{24}\\ -{\mathrm{\&alpha;}}_{13}& -{\mathrm{\&alpha;}}_{14}& {\mathrm{\&alpha;}}_{11}& -{\mathrm{\&alpha;}}_{12}\\ {\mathrm{\&alpha;}}_{23}& {\mathrm{\&alpha;}}_{24}& -{\mathrm{\&alpha;}}_{12}& {\mathrm{\&alpha;}}_{22}\end{array}\right]\end{array}---\left(18\right)$

进一步简化式(18)，得到的优化表达式

$\begin{array}{c}\tilde{B}=2\left[\begin{array}{cccc}{\mathrm{\&alpha;}}_{11}& {\mathrm{\&alpha;}}_{12}& 0& 0\\ {\mathrm{\&alpha;}}_{12}& {\mathrm{\&alpha;}}_{22}& 0& 0\\ 0& 0& {\mathrm{\&alpha;}}_{33}& -{\mathrm{\&alpha;}}_{34}\\ 0& 0& -{\mathrm{\&alpha;}}_{34}& {\mathrm{\&alpha;}}_{44}\end{array}\right]+2\left[\begin{array}{cccc}{\mathrm{\&alpha;}}_{33}& {\mathrm{\&alpha;}}_{34}& 0& 0\\ {\mathrm{\&alpha;}}_{34}& {\mathrm{\&alpha;}}_{44}& 0& 0\\ 0& 0& {\mathrm{\&alpha;}}_{11}& -{\mathrm{\&alpha;}}_{12}\\ 0& 0& -{\mathrm{\&alpha;}}_{12}& {\mathrm{\&alpha;}}_{22}\end{array}\right]\\ =2\left[\begin{array}{cccc}{\mathrm{\&alpha;}}_{11}+{\mathrm{\&alpha;}}_{33}& {\mathrm{\&alpha;}}_{12}+{\mathrm{\&alpha;}}_{34}& 0& 0\\ {\mathrm{\&alpha;}}_{12}+{\mathrm{\&alpha;}}_{34}& {\mathrm{\&alpha;}}_{22}+{\mathrm{\&alpha;}}_{44}& 0& 0\\ 0& 0& {\mathrm{\&alpha;}}_{11}+{\mathrm{\&alpha;}}_{33}& -{\mathrm{\&alpha;}}_{12}-{\mathrm{\&alpha;}}_{34}\\ 0& 0& -{\mathrm{\&alpha;}}_{12}-{\mathrm{\&alpha;}}_{34}& {\mathrm{\&alpha;}}_{22}+{\mathrm{\&alpha;}}_{44}\end{array}\right]\end{array}---\left(19\right)$

由式(19)看出，只需要求出α₁₁、α₂₂、α₃₃、α₄₄、α₁₂和α₃₄，即得到

同样，得到式(14)的向量表达式

令其中，c₁、c₂、c₃、c₄分别表示各个维度的值，那么由式(12)得到令其中，w₁表示w_ij的实部，w₂表示元素w_ji的实部，w₃表示w_ij的虚部，w₄表示w_ji的虚部，得到：

$2\left[\begin{array}{cccc}{\mathrm{\&alpha;}}_{11}+{\mathrm{\&alpha;}}_{33}& {\mathrm{\&alpha;}}_{12}+{\mathrm{\&alpha;}}_{34}& 0& 0\\ {\mathrm{\&alpha;}}_{12}+{\mathrm{\&alpha;}}_{34}& {\mathrm{\&alpha;}}_{22}+{\mathrm{\&alpha;}}_{44}& 0& 0\\ 0& 0& {\mathrm{\&alpha;}}_{11}+{\mathrm{\&alpha;}}_{33}& -{\mathrm{\&alpha;}}_{12}-{\mathrm{\&alpha;}}_{34}\\ 0& 0& -{\mathrm{\&alpha;}}_{12}-{\mathrm{\&alpha;}}_{34}& {\mathrm{\&alpha;}}_{22}+{\mathrm{\&alpha;}}_{44}\end{array}\right]\left[\begin{array}{c}{w}_1\\ w_2\\ w_3\\ w_4\end{array}\right]=\left[\begin{array}{c}{c}_1\\ c_3\\ c_2\\ c_4\end{array}\right]---\left(21\right)$

然后利用i∈{1,…,N-1}、j∈{i+1,…,N}时对应的所有值，计算得到第n次迭代后包含的代价函数的绝对值；设定收敛条件为：

$\left|\min L\left(\widehat{w_{ij}}\right)\right|\&le;\mathrm{\&epsiv;}$

其中，ε表示设定的趋于无穷小的极小数，本发明设定ε等于10^-6；将所述第n次迭代后包含的代价函数的绝对值与ε进行比较；若所述第n次迭代后包含的代价函数的绝对值大于ε，则算法未收敛，令n加1，返回步骤3；

若所述第n次迭代后包含的代价函数的绝对值小于等于ε，则算法达到收敛，此时依次计算得到第n次迭代后第k个雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵W_(n)、第k个雷达目标矩阵经过第n次迭代后得到的N×N维解混迭矩阵V_(n+1)和第k个雷达目标矩阵经过第n次迭代后的雷达目标矩阵对角化结构具体为：

利用i∈{1,…,N-1}、j∈{i+1,…,N}时对应的所有值，并根据w_ij＝w₁+j*w₃、w_ji＝w₂+j*w₄计算得到w_ij和w_ji各自对应的所有值，进而得到第n次迭代后第k个雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵W_(n)。

根据V_(n+1)＝(I+W_(n))V_(n)计算得到第k个雷达目标矩阵经过第_n次迭代后得到的N×N维解混迭矩阵V_(n+1)；根据计算得到第k个雷达目标矩阵经过第n次迭代后的雷达目标矩阵对角化结构并将所述第k个雷达目标矩阵经过第n次迭代后的雷达目标矩阵对角化结构作为第k个雷达目标矩阵的联合对角化结构；然后令k加1，同时初始化n为1，返回步骤3；直到k＝K，进而得到K个雷达目标矩阵各自对应的联合对角化结构，用于解决K个雷达目标矩阵中的盲源分离问题。

其中，表示w_ij的实部，表示w_ji的实部，表示w_ij的虚部，表示w_ji的虚部，w_ij表示W_(n)的第i行、第j列元素，w_ji表示W_(n)的第j行、第i列元素，W_(n)表示第n次迭代后第k个雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵，上标T表示转置。

下面通过仿真实验对本发明效果做进一步验证。

首先给出两个性能指标，其一为矩阵VC^kV^H中所有非对角线元素的平方和，k∈{1,…,K}，并命名对角化误差(Diagonalization Error,DE)为：

$DE=\underset{k=1}{\overset{K}{\&Sigma;}}\underset{i=1}{\overset{N}{\&Sigma;}}\underset{j=1,j\&NotEqual;i}{\overset{N}{\&Sigma;}}|{\left({\mathrm{VC}}^k{V}^H\right)}_{ij}{|}^2$

用于反映CVFFDIAG算法所得代价函数的变化趋势；其二为全局拒绝水平(GRL)参数；在盲源分离领域中，该全局拒绝水平(GRL)参数被公认为一个科学且重要的性能参数，用来表达通过估计得到的解混迭矩阵与混迭矩阵之间的差别。

(一)实验条件

实验一：构造无噪声雷达目标矩阵组来验证本发明方法的收敛性能；

实验二：在雷达目标矩阵组存在噪声的条件下验证本发明方法的收敛性能；

实验三：合成带噪雷达目标矩阵组，对比本发明方法和ACDC算法、FAJD算法的性能差异。

(二)实验内容

实验一：给出N×N维目标方阵组C^k＝AΛ^kA^H，其中K＝15，N＝5，随机产生复值矩阵分别为混迭矩阵A和第k个雷达目标矩阵的对角阵Λ^k，k∈{1,…,K}；给定N×N维解混迭矩阵V的初值：V₍₁₎＝I。运行20次后，仿真结果如图2A和图2B所示；图2A是未加入噪声的全局拒绝水平随迭代次数变化的曲线示意图，图2B是未加入噪声的对角化误差随迭代次数变化的曲线示意图。

实验二：给出含噪声的第k个N×N维雷达目标矩阵C^k＝AΛ^kA^H+ΔC^k，k∈{1,…,K}；其中K＝15，N＝5；第k个N×N维雷达目标矩阵的误差矩阵ΔC^k为复值矩阵，为了表征噪声矩阵的强度，将含噪声的第k个N×N维雷达目标矩阵C^k中不含噪声矩阵AΛ^kA^H和噪声矩阵ΔC^k的比表示为NER：

$NER=10{\log }_{10}\frac{\left|\right|{\mathrm{A\&Lambda;}}^k{A}^H|{|}_F^2}{\left|\right|{\mathrm{\&Delta;C}}^k|{|}_F^2}$

并产生分别满足NER＝10dB、15dB、20dB的噪声矩阵ΔC^k，将本发明方法分别在三种不同的NER条件下独立实验20次后，仿真结果如图3A和图3B所示，图3A是加入噪声的全局拒绝水平随迭代次数变化的曲线示意图，图3B是加入噪声的平均全局拒绝水平随迭代次数变化的曲线示意图；其中，图3A和图3B中的点线：NER＝10dB，虚线：NER＝15dB，实线：NER＝20dB。

实验三：给出第k个N×N维雷达目标矩阵C^k＝AΛ^kA^H+ΔC^k，其中K＝15，N＝5；本实验中，由于ACDC算法中要求雷达目标矩阵组必须满足Hermitian的要求，所以对角阵Λ^k(k＝1,…,K)均为实值矩阵，且第k个N×N维雷达目标矩阵中的噪声矩阵ΔC^k同样也满足Hermitian要求，产生分别满足NER＝10dB、15dB、20dB噪声矩阵ΔC^k；将本发明方法分别在三种不同的NER条件下独立实验20次后，仿真结果如图4所示，图4是分别使用DC算法、FAJD算法和本发明方法得到的全局拒绝水平随迭代次数变化的曲线示意图；其中，点线：ACDC算法，虚线：FAJD算法，实线：本发明方法。

(三)结果分析

实验一：图2A和图2B表明，在无噪声的条件下，本发明方法能够迅速收敛，而且其收敛值与实际值的差别很小。

实验二：图3A和3B表明，在带噪的条件下，本发明方法在10次迭代内便收敛到一个很小的值，说明本发明方法即使在带噪的条件下，也同样能够迅速收敛且具有很好地收敛性能。

实验三：图4的比较结果显示，本发明方法的收敛速度最快，FAJD算法次之，ACDC算法的收敛速度最慢；而从收敛后的趋势来看，本发明方法要略差于ACDC算法、然而优于FAJD算法，纵然在噪声较大的情况下，如NER＝10dB，本发明方法的收敛误差也低于-27dB，而且ACDC算法要求雷达目标矩阵组必须满足Hermitian要求，所以总的来说，本发明方法优于ACDC算法和FAJD算法，且具有快速地收敛速度和很好地收敛特性。

综上所述，仿真实验验证了本发明算法的正确性，有效性和可靠性。

显然，本领域的技术人员可以对本发明进行各种改动和变型而不脱离本发明的精神和范围；这样，倘若本发明的这些修改和变型属于本发明权利要求及其等同技术的范围之内，则本发明也意图包含这些改动和变型在内。

Claims (7)

1.一种基于FFDIAG算法的雷达目标盲源分离方法，其特征在于，包括以下步骤：

步骤1，分别设定雷达目标矩阵总个数K以及N×N维复矩阵集合C^N×N，并计算得到基于F-范数的雷达目标矩阵代价函数L(V)；其中，V表示设定的解混迭矩阵，K和N分别为自然数；

步骤2，初始化：令n为迭代次数，且初始值为1；k表示第k个雷达目标矩阵，且k的初始值为1，V₍₁₎表示N×N维单位阵，W₍₀₎表示N×N维零矩阵；

步骤3，计算得到第k个雷达目标矩阵经过第n次迭代后得到的N×N维解混迭矩阵V_(n+1)，并将V_(n+1)写成与W_(n)有关的代价函数形式L(W_(n))；其中，W_(n)表示第n次迭代后第k个雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵，N为自然数；

步骤4，根据与W_(n)有关的代价函数形式L(W_(n))，计算得到关于的代价函数其中，令 表示w_ij的实部，表示w_ji的实部，表示w_ij的虚部，表示w_ji的虚部，w_ij表示W_(n)的第i行、第j列元素，w_ji表示W_(n)的第j行、第i列元素，W_(n)表示第n次迭代后第k个雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵，上标T表示转置；

步骤5，将关于的代价函数对求导，并令导数等于零，得到的表达式然后令i的初值为1，j∈{i+1,…,N}，得到的N-1个对应值；再令i加1，j∈{i+1,…,N}，得到的N-2个对应值；直到i＝N-1，j∈{i+1,…,N}，得到的1个对应值；

利用i∈{1,…,N-1}、j∈{i+1,…,N}时对应的所有值，计算得到第n次迭代后包含的代价函数的绝对值并设定收敛条件：

$\left|\min L\left({\hat{w}}_{ij}\right)\right|\&le;\mathrm{\&epsiv;}$

其中，ε表示设定的趋于无穷小的极小数；将所述第n次迭代后包含的代价函数的绝对值与ε进行比较；若所述第n次迭代后包含的代价函数的绝对值大于ε，则令n加1，返回步骤3；

若所述第n次迭代后包含的代价函数的绝对值小于等于ε，则依次计算得到第n次迭代后第k个雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵W_(n)、第k个雷达目标矩阵经过第n次迭代后得到的N×N维解混迭矩阵V_(n+1)和第k个雷达目标矩阵经过第n次迭代后的雷达目标矩阵对角化结构并将所述第k个雷达目标矩阵经过第n次迭代后的雷达目标矩阵对角化结构作为第k个雷达目标矩阵的联合对角化结构；然后令k加1，同时初始化n为1，返回步骤3；直到k＝K，进而得到K个雷达目标矩阵各自对应的联合对角化结构，用于解决K个雷达目标矩阵中的盲源分离问题。

2.如权利要求1所述的一种基于FFDIAG算法的雷达目标盲源分离方法，其特征在于，在步骤1中，所述基于F-范数的雷达目标矩阵代价函数L(V)，其表达式为：

$L\left(V\right)=\underset{V\&Element;C^{N\&times;N}}{min}\underset{k=1}{\overset{K}{\&Sigma;}}off\left({\mathrm{VC}}^k{V}^H\right)$

其中，min表示求取最小值操作，C^N×N表示N×N维复矩阵集合，k∈{1,…,K}，K表示雷达目标矩阵总个数，V表示设定的解混迭矩阵，off(·)表示矩阵中所有非对角线元素的F-范数之和，C^k表示第k个雷达目标矩阵，上标H表示共轭转置，K和N分别为自然数。

3.如权利要求1所述的一种基于FFDIAG算法的雷达目标盲源分离方法，其特征在于，在步骤3中，所述第k个雷达目标矩阵经过第n次迭代后得到的N×N维解混迭矩阵V_(n+1)和所述与W_(n)有关的代价函数形式L(W_(n))，其表达式为：

V_(n+1)＝(I+W_(n))V_(n)

其中，I表示N×N维单位阵，V_(n)表示第k个雷达目标矩阵经过第n-1次迭代后得到的N×N维解混迭矩阵，min表示求取最小值操作，C^N×N表示N×N维复矩阵集合，W_(n)表示第n次迭代后第k个雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵，表示的对角部分，表示的非对角部分，表示第k个雷达目标矩阵经过第n-1次迭代后的雷达目标矩阵，上标H表示共轭转置。

4.如权利要求1所述的一种基于FFDIAG算法的雷达目标盲源分离方法，其特征在于，在步骤4中，所述关于的代价函数其表达式为：

其中，C^N×N表示N×N维复矩阵集合，D^k表示第k个雷达目标矩阵C^k的对角部分，E^k表示第k个雷达目标矩阵C^k的非对角部分，w_ij表示W_(n)的第i行、第j列元素，表示E^k的第i行、第j列元素，表示D^k对角线上的第i个元素，w_ji表示W_(n)的第j行、第i列元素，W_(n)表示第n次迭代后第k个雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵；表示E^k的第j行、第i列元素，表示D^k对角线上的第j个元素，D^k表示第k个雷达目标矩阵C^k的对角部分，E^k表示第k个雷达目标矩阵C^k的非对角部分。

5.如权利要求1所述的一种基于FFDIAG算法的雷达目标盲源分离方法，其特征在于，在步骤5中，所述将关于的代价函数对求导，并令导数等于零，还包括：

将关于的代价函数对进行求导，并令其导数等于零，得到：

其中， 上标T表示转置，表示w_ij的实部，表示w_ji的实部，表示w_ij的虚部，表示w_ji的虚部，w_ij表示W_(n)的第i行、第j列元素，w_ji表示W_(n)的第j行、第i列元素，W_(n)表示第n次迭代后第k个雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵；表示的实部，表示的实部，表示的虚部，表示的虚部，表示D^k对角线上的第i个元素，表示D^k对角线上的第j个元素，D^k表示第k个雷达目标矩阵C^k的对角部分。

6.如权利要求1所述的一种基于FFDIAG算法的雷达目标盲源分离方法，其特征在于，在步骤5中，所述的表达式为：

其中， 上标T表示转置，表示w_ij的实部，表示w_ji的实部，表示w_ij的虚部，表示w_ji的虚部，w_ij表示W_(n)的第i行、第j列元素，w_ji表示W_(n)的第j行、第i列元素，W_(n)表示第n次迭代后第k个雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵；表示的实部，表示的实部，表示的虚部，表示的虚部，表示D^k对角线上的第i个元素，表示D^k对角线上的第j个元素，D^k表示第k个雷达目标矩阵C^k的对角部分，K表示雷达目标矩阵总个数，表示的实部，表示的虚部，表示E^k的第j行、第i列元素，E^k表示第k个雷达目标矩阵C^k的非对角部分。

7.如权利要求1所述的一种基于FFDIAG算法的雷达目标盲源分离方法，其特征在于，在步骤5中，所述依次计算得到第n次迭代后第k个雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵W_(n)、第k个雷达目标矩阵经过第n次迭代后得到的N×N维解混迭矩阵V_(n+1)和第k个雷达目标矩阵经过第n次迭代后的雷达目标矩阵对角化结构其过程为：

利用i∈{1,…,N-1}、j∈{i+1,…,N}时对应的所有值，并根据w_ij＝w₁+j*w₃、w_ji＝w₂+j*w₄计算得到w_ij和w_ji各自对应的所有值，进而得到第n次迭代后第k个雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵W_(n)；根据V_(n+1)＝(I+W_(n))V_(n)计算得到第k个雷达目标矩阵经过第n次迭代后得到的N×N维解混迭矩阵V_(n+1)；根据计算得到第k个雷达目标矩阵经过第n次迭代后的雷达目标矩阵对角化结构

其中，w_ij表示W_(n)的第i行、第j列元素，w_ji表示W_(n)的第j行、第i列元素，w₁表示w_ij的实部，w₂表示元素w_ji的实部，w₃表示w_ij的虚部，w₄表示w_ji的虚部，W_(n)表示第n次迭代后第k个雷达目标矩阵中主对角线元素被限制为零的N×N维矩阵。