CN103163505A - 基于jade的时变窄带干扰抑制方法 - Google Patents

Abstract

本发明公开了一种基于JADE的时变窄带干扰抑制方法，其步骤包括：(1)雷达录取回波；(2)判断是否存在窄带干扰；(3)预处理；(4)分离基信号；(5)提取干扰信号；(6)剔除干扰信号；(7)判断是否遍历完所有方位回波；(8)结束遍历。本发明通过采用近似联合对角化JADE方法直接在信号数据域提取窄带干扰，充分利用雷达信号与干扰信号统计特性的差异性，能够有效地抑制时变的窄带干扰，保留有用回波信号，克服了小波变换法耗时多、复数经验模态分解法对干扰信号分离不彻底的不足，具有信号损失小、抑制效果好的优点。

Description

基于JADE的时变窄带干扰抑制方法

技术领域

本发明属于信号处理技术领域，更进一步涉及雷达信号处理技术领域中的一种基于近似联合对角化(Joint Approximate Diagonalization of Eigen-matrices，JADE)的时变窄带干扰抑制方法。本发明可在有效地保留有用雷达回波信号的同时，尽可能地抑制针对雷达的时变窄带干扰，从而实现对雷达回波信号清晰成像。

背景技术

窄带干扰是针对合成孔径雷达系统的主要干扰方式之一。所谓窄带干扰是指干扰信号带宽与合成孔径雷达发射信号带宽之比很小(比如1％)的信号。窄带干扰可分为非人为干扰和人为干扰两类。非人为干扰主要包括与SAR处于同一频段的通讯设备、电视网和其它低频波段辐射设备。人为窄带干扰具有干扰功率大，针对性强，干扰实现简单，形式多变等特点，最常见的有强脉冲干扰和频率引导“瞄准式”干扰。这些干扰的存在降低了SAR回波数据的信干噪比，尤其是在干扰功率较大的情况下，会使得图像中出现亮线并且使图像变模糊，从而严重降低了SAR图像质量。因此，如何对窄带干扰形式进行有效的识别和抑制是目前SAR干扰处理技术的重要课题。

张双喜等人在文献“利用小波变换进行合成孔径雷达窄带干扰抑制算法”(西安电子科技大学学报自然科学版，vol.38，no.5，pp.108-114，2011年)中提出一种基于小波变换的窄带干扰抑制方法。该方法首先将回波信号的频谱变换到小波域中，再把窄带干扰所对应的小波系数检测出来并通过设计具体的时频域滤波器去除干扰，从而实现窄带干扰抑制。但是该方法存在的不足之处是，在不同的尺度下需要对选定提取的窄带干扰的阈值进行恒虚警率检测，从而对小波变换中的小波系数检测需要消耗较多的时间。

西安电子科技大学在其申请的专利“基于复数经验模态分解的时变窄带干扰抑制方法”(专利申请号201210001678.2，公开号CN102520396A)中公开了一种利用复数经验模态分解的时变窄带干扰抑制方法。该方法将信号分解成一系列本征模态函数，在抑制窄带干扰的同时充分保留了信号本身的非线性和非平稳特征，但是该方法仍存在的不足是，对局部信号分析处理未考虑干扰信号的全局统计特性，从而导致干扰信号分离不彻底。

发明内容

本发明的目的在于克服现有技术的不足，提出一种基于JADE的时变窄带干扰抑制方法。该方法弥补了小波变换法耗时比较多、复数经验模态分解法对干扰信号分离不彻底的不足，充分利用雷达信号与干扰信号统计特性的差异性，根据JADE分解方法直接在信号数据域提取窄带干扰并进行抑制，能够保留有用回波信号，有效地抑制时变的窄带干扰。

实现本发明的基本思路是：首先，在频域对回波数据进行窄带干扰定性识别，对存在干扰的回波进行时间平滑，并对信号进行白化处理；其次，利用SAR回波和窄带干扰信号的独立性，对回波进行JADE分解，得到一系列的基信号，通过能量阈值检测实现干扰基信号识别；最后，采用最小二乘方法获得对干扰基信号的复包络的精确估计后，重构干扰信号，并从原始回波中消去干扰信号，实现方位时变窄带干扰抑制。

本发明的具体步骤如下：

(1)雷达录取回波，得到以距离为行向量，以方位为列向量的回波数据；

(2)判断回波数据是否存在窄带干扰：

2a)对回波数据的列向量进行傅里叶变换，得到一维距离像；

2b)按照下式计算一维距离像的相对幅度比：

$R=\frac{\max \left(\right|P\left|\right)}{\mathrm{mean}\left(\right|P\left|\right)}$

其中，R表示一维距离像的相对幅度比，max(·)表示求最大值运算符号，|·|表示求模值运算符号，P表示一维距离像，mean(·)表示求均值运算符号；

2c)判断一维距离像的相对幅度比是否大于或等于幅度阈值，若是，则执行步骤(3)；否则，执行步骤(7)；

(3)预处理：

3a)以定长在回波数据的列向量中滑动选取数据，以每次所选取的数据为行向量，组成数据矩阵；

3b)分别计算数据矩阵每行的均值，将数据矩阵的各行数据分别与其所在行对应的均值相减，并将差值依次作为行向量中的元素，组成零均值数据矩阵；

3c)计算零均值数据矩阵的协方差矩阵，并对协方差矩阵进行特征值分解，得到特征值对角矩阵和特征向量矩阵；

3d)按照下式计算白化矩阵：

W＝Λ^-1/2Q^H

其中，W表示白化矩阵，Λ表示特征值对角矩阵，Λ^-1/2表示对Λ进行开方求逆得到的矩阵，Q表示特征向量矩阵，Q^H表示Q的共轭转置矩阵；

3e)将白化矩阵与零均值数据矩阵相乘，得到白化观测矩阵；

(4)分离基信号：

4a)按照下式计算白化观测矩阵的四阶协方差矩阵：

$C=\frac{1}{N}\underset{k=1}{\overset{N}{\mathrm{\Σ}}}\left[Z_k^H{Z}_k{Z}_k{Z}_k^H\right]-2\mathrm{RR}-\mathrm{tr}\left(R\right)R$

其中，C表示白化观测矩阵的四阶协方差矩阵，N表示白化观测矩阵的列数，

表示对第1项到第N项的求和运算符号，k表示白化观测矩阵列向量的编号，Z_k表示白化观测矩阵第k列向量，表示Z_k的共轭转置向量，R表示白化观测矩阵的自相关矩阵，tr(R)表示对自相关矩阵R求矩阵迹得到的数值；

4b)对四阶协方差矩阵进行特征值分解，得到一组特征向量，依次将每一个特征向量与自身的共轭转置向量相乘，得到各特征向量对应的转换矩阵；

4c)按照下式计算四阶累积量矩阵：

$C_p=\frac{1}{N}\underset{k=1}{\overset{N}{\mathrm{\Σ}}}\left[Z_k^H{E}_p{Z}_k{Z}_k{Z}_k^H\right]-{\mathrm{RE}}_{p}R-\mathrm{tr}\left(E_{p}R\right)R-{\mathrm{RE}}_p^{H}R$

其中，C_P表示第p个特征向量对应的四阶累积量矩阵，p表示特征向量的标号，p的取值范围为[1，L]，L表示步骤3a)中选取的定长，N表示白化观测矩阵的列数，

表示对第1项到第N项的求和运算符号，k表示白化观测矩阵列向量的编号，Z_k表示白化观测矩阵第k列向量，表示Z_k的共轭转置向量，E_p表示第p个特征向量对应的转换矩阵，R表示白化观测矩阵的自相关矩阵，tr(·)表示求矩阵迹运算符号，

表示E_p的共轭转置矩阵；

4d)按照联合对角化法计算Givens旋转矩阵；

4e)利用JADE原理计算基信号矩阵；

(5)提取干扰信号：

5a)依次将基信号矩阵中的行向量作为基信号，按照下式计算各基信号的峰值能量比：

${\mathrm{\τ}}_m=\frac{\max \left[\mathrm{FT}\left(y_m\right)\right]}{E\left[\mathrm{FT}\left(y_m\right)\right]}$

其中，τ_m表示第m个基信号的峰值能量比，m表示基信号矩阵的行数，max(·)表示计算最大值运算符号，FT(·)表示傅里叶变换运算符号，y_m表示基信号矩阵的第m行向量，E(·)表示求均值运算符号；

5b)将峰值能量比大于或等于能量阈值的基信号作为干扰基信号；

5c)以干扰基信号为行向量，组成干扰信号矩阵；

(6)剔除干扰信号：

6a)按照最小二乘法估计干扰信号的复包络，得到复包络向量；

6b)将复包络向量和步骤5c)得到的干扰信号矩阵相乘，将乘积作为重构的干扰信号；

6c)将回波数据的列向量与重构的干扰信号相减，得到剔除干扰后的回波数据的列向量；

(7)判断是否遍历完所有回波数据的列向量，若遍历完，则执行步骤(8)；否则，执行步骤(2)；

(8)获取剔除干扰后的回波数据。

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

第一，本发明通过在频域对窄带干扰进行定性识别，不需要估计任何参数，并对数据进行白化处理，降低数据矩阵的维数，克服了现有技术中小波变换方法需要精确估计干扰频率，多次设置阈值导致运算量大的不足，使得本发明具有运算量较少，效率较高的优点。

第二，本发明通过采用JADE分解方法，充分利用回波信号和干扰信号的统计特性差异，按照独立性度量指标将干扰信号与回波信号独立的分离出来，克服了现有技术中基于复数经验模态分解方法仅考虑局部统计特性，而造成干扰抑制不充分的缺点，使得本发明具有信号损失小，抑制效果好的优点。

附图说明

图1为本发明的流程图；

图2为本发明的仿真图。

具体实施方式

下面结合附图1，对本发明具体实施方式作进一步的详细描述。

步骤1，获取场景的SAR回波，雷达以脉冲重复频率发射并接收脉冲，得到以距离为行向量，以方位为列向量的SAR回波数据。

步骤2，判断回波数据是否存在窄带干扰。

对回波数据的列向量进行傅里叶变换，得到一维距离像；按照下式计算一维距离像的相对幅度比：

$R=\frac{\max \left(\right|P\left|\right)}{\mathrm{mean}\left(\right|P\left|\right)}$

其中，R表示一维距离像的相对幅度比，max(·)表示求最大值运算符号，|·|表示求模值运算符号，P表示一维距离像，mean(·)表示求均值运算符号；判断一维距离像的相对幅度比是否大于或等于幅度阈值，若是，则执行步骤3；否则，执行步骤7。

步骤3，回波数据预处理。

3a)以定长在回波数据的列向量中滑动选取数据，以每次所选取的数据为行向量，组成数据矩阵，其中，定长的限定如下：

$n\≤L\≤\left(\frac{N-1}{2}\right)$

其中，n表示用户选定的窄带干扰分量个数，通常取1～32之间，L表示定长，N表示距离维采样点数目。

3b)分别计算数据矩阵每行的均值，将数据矩阵的各行数据分别与其所在行对应的均值相减，并将差值依次作为行向量中的元素，组成零均值数据矩阵。

3c)计算零均值数据矩阵的协方差矩阵，并对协方差矩阵进行特征值分解，得到特征值对角矩阵和特征向量矩阵，其中，特征值分解方法如下：

C_X＝QΛQ^H

其中，C_X表示零均值数据矩阵的协方差矩阵，Q表示协方差矩阵的特征向量矩阵，Λ表示协方差矩阵的特征值对角矩阵，Q^H表示Q的共轭转置矩阵。

3d)按照下式计算白化矩阵：

W＝Λ^-1/2Q^H

其中，W表示白化矩阵，Λ表示特征值对角矩阵，Λ_-1/2表示对Λ进行开方求逆得到的矩阵，Q表示特征向量矩阵，Q^H表示Q的共轭转置矩阵。

3e)将白化矩阵与步骤3b)得到的零均值数据矩阵相乘，得到白化观测矩阵。

步骤4，利用JADE分离基信号。

4a)按照下式计算白化观测矩阵的四阶协方差矩阵：

$C=\frac{1}{N}\underset{k=1}{\overset{N}{\mathrm{\Σ}}}\left[Z_k^H{Z}_k{Z}_k{Z}_k^H\right]-2\mathrm{RR}-\mathrm{tr}\left(R\right)R$

其中，C表示白化观测矩阵的四阶协方差矩阵，N表示白化观测矩阵的列数，

表示对第1项到第N项的求和运算符号，k表示白化观测矩阵列向量的编号，Z_k表示白化观测矩阵第k列向量，表示Z_k的共轭转置向量，R表示白化观测矩阵的自相关矩阵，tr(R)表示对自相关矩阵R求矩阵迹得到的数值。

4b)对四阶协方差矩阵按照下式进行特征值分解，得到一组特征向量，依次将每一个特征向量与自身的共轭转置向量相乘，得到各特征向量对应的转换矩阵：

C＝UDU^H

其中，C表示步骤4a)得到的四阶协方差矩阵，U表示四阶协方差矩阵的特征向量矩阵，U＝[u₁，u₂，...，u_L]，u₁，u₂，...，u_L表示特征向量，L表示步骤3a)中限定的定长，D表示四阶协方差矩阵的特征值对角矩阵，D＝diag[λ₁，λ₂，...，λ_L]，diag[λ₁，λ₂，...，λ_L]表示以特征值λ₁，λ₂，...，λ_L为对角线元素的对角矩阵，U^H表示U的共轭转置矩阵。

4c)按照下式计算四阶累积量矩阵：

$C_p=\frac{1}{N}\underset{k=1}{\overset{N}{\mathrm{\Σ}}}\left[Z_k^H{E}_p{Z}_k{Z}_k{Z}_k^H\right]-{\mathrm{RE}}_{p}R-\mathrm{tr}\left(E_{p}R\right)R-{\mathrm{RE}}_p^{H}R$

其中，C_P表示第p个特征向量对应的四阶累积量矩阵，p表示特征向量的标号，p的取值范围为[1，L]，L表示步骤3a)中选取的定长，N表示白化观测矩阵的列数，

表示对第1项到第N项的求和运算符号，k表示白化观测矩阵列向量的编号，Z_k表示白化观测矩阵第k列向量，

表示Z_k的共轭转置向量，E_p表示第p个特征向量对应的转换矩阵，R表示白化观测矩阵的自相关矩阵，tr(·)表示求矩阵迹运算符号，

表示E_p的共轭转置矩阵。

4d)按照联合对角化法计算Givens旋转矩阵，其中，联合对角化法如下：

$G=\arg \min \left(\underset{p=1}{\overset{L}{\mathrm{\Σ}}}\mathrm{off}\right[G^H{C}_{P}G\left]\right)$

其中，G表示Givens旋转矩阵，argmin(·)表示计算使目标函数取最小值时对应变量值的运算符号，

表示对第1项到第L项求和的运算符号，p表示步骤4b)得到的特征向量的序号，L表示步骤3a)中选取的定长，off(·)表示对矩阵非对角线元素求和运算符号，G^H表示Givens旋转矩阵G的共轭转置矩阵，C_P表示步骤4c)得到的第p个特征向量对应的四阶累积量矩阵。

4e)利用JADE原理计算基信号矩阵，其中，JADE原理如下：

S＝G^HWX

其中，S表示基信号矩阵，G表示步骤4d)得到的Givens旋转矩阵，G^H表示G的共轭转置矩阵，W表示步骤3d)得到的白化矩阵，X表示步骤3e)得到的白化观测矩阵。

步骤5，提取干扰信号。

依次将基信号矩阵中的行向量作为基信号，按照下式计算各基信号的峰值能量比：

${\mathrm{\τ}}_m=\frac{\max \left[\mathrm{FT}\left(y_m\right)\right]}{E\left[\mathrm{FT}\left(y_m\right)\right]}$

其中，τ_m表示第m个基信号的峰值能量比，m表示基信号矩阵的行数，max(·)表示计算最大值运算符号，FT(·)表示傅里叶变换运算符号，y_m表示基信号矩阵的第m行向量，E(·)表示求均值运算符号；将峰值能量比大于或等于能量阈值的基信号作为干扰基信号，其中，能量阈值限定在2～10之间；以干扰基信号为行向量，组成干扰信号矩阵；

步骤6，剔除干扰信号。

6a)按照最小二乘法估计干扰信号的复包络，得到复包络向量，其中，最小二乘法如下：

A＝(I^HI)^-1I^Hd

其中，A表示干扰信号的复包络向量，I表示干扰信号矩阵，I^H表示干扰信号矩阵I的共轭转置矩阵，(·)^-1表示矩阵求逆运算符号，d表示回波数据的列向量。

6b)将复包络向量和步骤5得到的干扰信号矩阵相乘，将乘积作为重构的干扰信号，其中重构的干扰信号可按下式进行计算：

J＝AI

其中，J表示重构的干扰信号，A表示步骤6a)得到的干扰信号的复包络向量，I表示步骤6a)得到的干扰信号矩阵。

6c)将回波数据的列向量与重构的干扰信号相减，得到剔除干扰后的回波数据的列向量，其中，剔除干扰后的回波数据的列向量可按下式进行计算：

F＝d-J

其中，F表示剔除干扰后的回波数据的列向量，d表示回波数据的列向量，J为步骤6b)得到的重构干扰信号。

步骤7，判断是否遍历完所有回波数据的列向量，若遍历完，则执行步骤8；否则，执行步骤2。

步骤8，获取剔除干扰后的回波数据。

下面结合附图2所示的仿真图对本发明的效果做进一步的描述。

1.仿真条件

附图2所示的仿真图是在MATLAB(R2010b)软件下进行的。由于雷达回波信号可以近似为线性调频信号，因此，仿真观测信号由一个线性调频信号和一个窄带干扰信号混合组成。原线性调频信号的调频率为6×10¹²Hz/s，干扰信号载频f_c为10MHz，采样频率为150MHz，观测时间为[-3.4，3.4]μs。

2.仿真内容

首先，对混合信号进行傅里叶变换，得到附图2(a)的距离维幅频特性图；然后，对混合信号进行WVD变换，得到附图2(b)的结果图；接着，采用JADE分解方法直接在信号数据域提取窄带干扰分量，并对窄带干扰分量进行WVD变换，得到附图2(c)的结果图；最后，重构干扰信号并从混合信号中剔除干扰信号，得到附图2(d)的结果图。

3.附图说明

附图2(a)是对混合信号进行傅里叶变换的距离维幅频特性图，其中，横坐标表示频域单元，纵坐标表示归一化幅度值，从图中可看出，由于窄带干扰的存在，线性调频信号频谱中存在明显的带内冲激干扰，在一个小频率段中幅值出现了突变，幅值比其它幅值高出了10倍左右。

附图2(b)是对混合信号进行WVD变换的结果图，其中，横坐标表示时间单元，纵坐标表示频率单元，由图中可看出，斜线表示线性调频信号的WVD变换结果，水平线表示干扰分量的WVD变换结果，图2(b)中的水平线与斜线交叉处存在阴影区域，这是由于观测信号为线性调频信号和干扰的混合，因此WVD变换结果中不仅存在自有项，还存在大量的交叉项，说明线性调频信号受到了干扰信号的影响。

附图2(c)是利用本文算法提出的干扰分量的WVD变换结果图，其中，横坐标表示时间单元，纵坐标表示频率单元。由图2(c)中可见，图2(c)中的水平线表示干扰分量的WVD变换结果，说明提取出来的干扰信号只有单频信息，时频信息很好的保留了下来。

附图2(d)是剔除干扰后的有用信号的WVD变换结果图，其中，横坐标表示时间单元，纵坐标表示频率单元。由图2(d)中可见，图2(d)中的斜线表示干扰分量的WVD变换结果，说明被恢复出来的信号时频域中已没有与窄带干扰的交叉项，只保留了原线性调频信号的时频信息。因此证明了本方法的有效性。

Claims (7)

1.基于JADE的时变窄带干扰抑制方法，包括如下步骤：

(1)雷达录取回波，得到以距离为行向量，以方位为列向量的回波数据；

(2)判断回波数据是否存在窄带干扰：

2a)对回波数据的列向量进行傅里叶变换，得到一维距离像；

2b)按照下式计算一维距离像的相对幅度比：

$R=\frac{\max \left(\right|\left|P\right|\left|\right)}{\mathrm{mean}\left(\right|\left|P\right|\left|\right)}$

其中，R表示一维距离像的相对幅度比，max(·)表示求最大值运算符号，||·||表示求模值运算符号，P表示一维距离像，mean(·)表示求均值运算符号；

2c)判断一维距离像的相对幅度比是否大于或等于幅度阈值，若是，则执行步骤(3)；否则，执行步骤(7)；

(3)预处理：

3a)以定长在回波数据的列向量中滑动选取数据，以每次所选取的数据为行向量，组成数据矩阵；

3b)分别计算数据矩阵每行的均值，将数据矩阵的各行数据分别与其所在行对应的均值相减，并将差值依次作为行向量中的元素，组成零均值数据矩阵；

3c)计算零均值数据矩阵的协方差矩阵，并对协方差矩阵进行特征值分解，得到特征值对角矩阵和特征向量矩阵；

3d)按照下式计算白化矩阵：

W＝Λ^-1/2Q^H

其中，W表示白化矩阵，Λ表示特征值对角矩阵，Λ^-1/2表示对Λ进行开方求逆得到的矩阵，Q表示特征向量矩阵，Q^H表示Q的共轭转置矩阵；

3e)将白化矩阵与步骤3b)得到的零均值数据矩阵相乘，得到白化观测矩阵；

(4)分离基信号：

4a)按照下式计算白化观测矩阵的四阶协方差矩阵：

$C=\frac{1}{N}\underset{k=1}{\overset{N}{\mathrm{\Σ}}}\left[Z_k^H{Z}_k{Z}_k{Z}_k^H\right]-2\mathrm{RR}-\mathrm{tr}\left(R\right)R$

其中，C表示白化观测矩阵的四阶协方差矩阵，N表示白化观测矩阵的列数，表示对第1项到第N项的求和运算符号，k表示白化观测矩阵列向量的编号，Z_k表示白化观测矩阵第k列向量，

表示Z_k的共轭转置向量，R表示白化观测矩阵的自相关矩阵，tr(R)表示对自相关矩阵R求矩阵迹得到的数值；

4b)对四阶协方差矩阵进行特征值分解，得到一组特征向量，依次将每一个特征向量与自身的共轭转置向量相乘，得到各特征向量对应的转换矩阵；

4c)按照下式计算四阶累积量矩阵：

$C_p=\frac{1}{N}\underset{k=1}{\overset{N}{\mathrm{\Σ}}}\left[Z_k^H{E}_p{Z}_k{Z}_k{Z}_k^H\right]-{\mathrm{RE}}_{p}R-\mathrm{tr}\left(E_{p}R\right)R-{\mathrm{RE}}_p^{H}R$

其中，C_P表示第p个特征向量对应的四阶累积量矩阵，p表示特征向量的标号，p的取值范围为[1，L]，L表示步骤3a)中选取的定长，N表示白化观测矩阵的列数，

表示对第1项到第N项的求和运算符号，k表示白化观测矩阵列向量的编号，Z_k表示白化观测矩阵第k列向量，

表示Z_k的共轭转置向量，E_p表示第p个特征向量对应的转换矩阵，R表示白化观测矩阵的自相关矩阵，tr(·)表示求矩阵迹运算符号，

表示E_p的共轭转置矩阵；

4d)按照联合对角化法计算Givens旋转矩阵；

4e)利用JADE计算基信号矩阵；

(5)提取干扰信号：

5a)依次将基信号矩阵中的行向量作为基信号，按照下式计算各基信号的峰值能量比：

${\mathrm{\τ}}_m=\frac{\max \left[\mathrm{FT}\left(y_m\right)\right]}{E\left[\mathrm{FT}\left(y_m\right)\right]}$

其中，τ_m表示第m个基信号的峰值能量比，m表示基信号矩阵的行数，max(·)表示计算最大值运算符号，FT(·)表示傅里叶变换运算符号，y_m表示基信号矩阵的第m行向量，E(·)表示求均值运算符号；

5b)将峰值能量比大于或等于能量阈值的基信号作为干扰基信号；

5c)以干扰基信号为行向量，组成干扰信号矩阵；

(6)剔除干扰信号：

6a)按照最小二乘法估计干扰信号的复包络，得到复包络向量；

6b)将复包络向量和步骤5c)得到的干扰信号矩阵相乘，将乘积作为重构的干扰信号；

6c)将回波数据的列向量与重构的干扰信号相减，得到剔除干扰后的回波数据的列向量；

(7)判断是否遍历完所有回波数据的列向量，若遍历完，则执行步骤(8)；否则，执行步骤(2)；

(8)获取剔除干扰后的回波数据。

2.根据权利要求1所述的基于JADE的时变窄带干扰抑制方法，其特征在于，步骤2c)中所述的幅度阈值限定在2～5之间。

3.根据权利要求1所述的基于JADE的时变窄带干扰抑制方法，其特征在于，步骤3a)中所述的定长限定如下：

$n\≤L\≤\left(\frac{N-1}{2}\right)$

其中，n表示用户选定的窄带干扰分量个数，通常取1～32之间，L表示定长，N表示距离维采样点数目。

4.根据权利要求1所述的基于JADE的时变窄带干扰抑制方法，其特征在于，步骤4d)中所述的联合对角化法如下：

$G=\arg \min \left(\underset{p=1}{\overset{L}{\mathrm{\Σ}}}\mathrm{off}\right[G^H{C}_{P}G\left]\right)$

其中，G表示Givens旋转矩阵，argmin(·)表示计算使目标函数取最小值时对应变量值的运算符号，表示对第1项到第L项求和的运算符号，p表示步骤4b)得到的特征向量的序号，L表示步骤3a)中选取的定长，off(·)表示对矩阵非对角线元素求和运算符号，G^H表示Givens旋转矩阵G的共轭转置矩阵，C_P表示步骤4c)得到的第p个特征向量对应的四阶累积量矩阵。

5.根据权利要求1所述的基于JADE的时变窄带干扰抑制方法，其特征在于，步骤4e)中所述的JADE原理如下：

S＝G^HWX

其中，S表示基信号矩阵，G表示步骤4d)得到的Givens旋转矩阵，G^H表示G的共轭转置矩阵，W表示步骤3d)得到的白化矩阵，X表示步骤3e)得到的白化观测矩阵。

6.根据权利要求1所述的基于JADE的时变窄带干扰抑制方法，其特征在于，步骤5b)中所述的能量阈值限定在2～10之间。

7.根据权利要求1所述的基于JADE的时变窄带干扰抑制方法，其特征在于，步骤6a)中所述的最小二乘法如下：

A＝(I^HI)^-1I^Hd

其中，A表示干扰信号的复包络向量，I表示干扰信号矩阵，I^H表示干扰信号矩阵I的共轭转置矩阵，(·)^-1表示矩阵求逆运算符号，d表示回波数据的列向量。

Images