CN103197325B - 一种基于变对角加载量的空时抗干扰方法 - Google Patents

Abstract

本发明公开了一种基于可变对角加载的空时抗干扰方法，包括以下几个步骤：步骤一，采集GPS卫星信号，信号经过采样后输入滤波器结构；步骤二，对信号进行相关相减多级维纳滤波器的前向分解，记录每级的接收信号和参考信号；步骤三，通过采样信号的协方差矩阵代替真实协方差矩阵；步骤四，通过矩阵求逆定理计算对角加载的范围，进而利用协方差对角阵元素关系求得需要加载的对角量；步骤五，加入对角量，进行相关相减多级维纳滤波器的后向综合，求出滤波器最优权矢量，对当前数据进行自适应滤波。本发明采用改进的相关相减多级维纳滤波器，相比多级维纳滤波器的运算量和存储量更低，适用于维数较大的系统。

Description

一种基于变对角加载量的空时抗干扰方法

技术领域

本发明提出一种基于变对角加载量的空时抗干扰方法，属于空时自适应阵列处理技术领域。

背景技术

随着导航卫星系统的广泛使用和信息化战争的发展，精确导航、制导成为了各方面的迫切需要。虽然GPS系统自身能够在某种程度上克服外界干扰，但是对于复杂的自然环境和特定的人为干扰，抗干扰技术的研究成为急需解决的关键问题。

自适应滤波是指利用前一时刻已获得的滤波器参数的结果，自动的调节当前时刻滤波器的参数，以适应信号和噪声随时间变化的统计特性，从而实现最优滤波。

自适应滤波主要有三个方面：时域自适应滤波，空域自适应滤波和频域自适应滤波。时域滤波是在时域内对信号进行处理，对窄带干扰、连续波干扰较为有效，能解决多径干扰和回波消除问题，但对宽带干扰滤除效果不佳。空域滤波主要有两种类型：零陷和波束形成。零陷技术对有用信号信息量的需求是最小的，它使信号的总输出功率最小化。波束形成技术最大化信号的输出信噪比。波束形成利用天线方向控制及各种天线阵列技术使波束成形达到最佳。天线方向控制用于减小对干扰的主波瓣、旁波瓣和后波瓣敏感性。空域滤波的不足之处主要有以下两个方面：一是如果阵元数为N，该阵最多能够产生的零陷数为N-1；二是如果一个干扰源离某个卫星的角度间隔很近，针对该干扰的空域零陷使该卫星的信号衰减得无法使用。频域滤波能够有效进行带通、带阻滤波，对抗窄带、连续波带内干扰，以及强带外干扰，但不能有效滤除宽带干扰和同一频率干扰。空时自适应抗干扰系统和空域自适应抗干扰系统具有相同的系统结构；不同之处是对每个阵元接收到的数据进行了L级延时，延迟单元相当于频域滤波器，可以在频域上对干扰进行抑制。它在不增加阵元的前提下，大大地增加了自由度，对于N个阵元，经过L次延迟的空时抗干扰系统，理论上能够抑制(N-1)×L个干扰，较之空域滤波抗干扰数目大大增加。

对角加载技术是目前广泛使用的改善波束形成的稳健方法，该技术可以在阵列响应失配的情况下，有效防止STAP主瓣发生畸变，改善自适应波束形成器的性能。该方法实现简单，减少了噪声特征值和扩散度，并且不增加运算量。但由于目前使用的加载方法多为固定值加载，在多变的工程应用中，加载量的大小很难确定，这样会使得抗干扰的性能在一定程度上下降，达不到改善滤波性能的目的。

发明内容

本发明的目的是为了提高空时自适应处理在小块拍下的稳健性，提出一种基于变对角加载量的空时抗干扰方法，本发明的可变对角加载量，首先利用矩阵求逆定理计算出加载量范围后，将加载量与阵元信号中的干扰和噪声进行相关，从而得到最佳的对角加载量。

本发明的一种基于变对角加载量的空时抗干扰方法，包括以下几个步骤：

步骤一：采集GPS卫星信号，信号经过采样后输入滤波器结构；

步骤二：对信号进行相关相减多级维纳滤波器的前向分解，记录每级的数据信号和参考信号；

步骤三：通过采样信号的协方差矩阵代替真实协方差矩阵；

步骤四：通过矩阵求逆定理计算对角加载的范围，进而利用干扰信号和噪声信号特征值的关系自适应的选择需要加载的对角量；

步骤五：加入对角量，进行相关相减多级维纳滤波器的后向综合，求出滤波器最优权矢量，对当前信号进行自适应滤波。

本发明的优点在于：

（1）本发明采用改进的相关相减多级维纳滤波器，相比多级维纳滤波器的运算量和存储量更低，适用于维数较大的系统；

（2）本发明使用变对角加载技术，既改善了小快拍下阵列响应失配现象，又避免了单一加载量不适合多变工作环境的情况，极大的提高了系统在小快拍下的稳健性。

附图说明

图1多级维纳滤波器结构图；

图2是相关相减多级维纳滤波器的结构图；

图3是本发明的方法流程图。

具体实施方式

下面将结合附图和实施例对本发明作进一步的详细说明。

图1是多级维纳滤波器结构，按照信号处理过程可分为前向分解和后向综合两个部分。前向分解是对采样后的卫星信号X(k)（图1中的X）进行正交投影，分解得到后一级的观测数据X₀(k)和参考信号d₀(k)，其中X(k)为N×1为向量，X₀(k)为(N-1)×1为向量，即观测数据维数降低1；以下的分解步骤如同第一步，每一级观测数据X_i(k)依次降低。h_i为第i级匹配滤波器，B_i为第i级阻塞矩阵，ε_i为第i级输出误差，其中i＝1,2,…,D≤N，N为阵元数。在最后一级分解中令X_D-1＝d_D＝ε_D。后向综合是由一组递推的标量综合得到维纳滤波器的输出误差信号ε_i，最后一级的观测数据与权系数w_D相乘并与参考信号做差值得到误差信号ε_D-1，依次迭代，最终得到系统输出误差值ε₀。

对多级维纳滤波器利用相关相减结构，得到图2的相关相减多级维纳滤波器，即CSA-MWF。相关相减多级维纳滤波器对阻塞矩阵进行特定要求，即通过等价结构，在第1级之后的分解过程中可以直接用如图2结构代替阻塞矩阵。图2中信号的处理过程与图1中的维纳滤波器相似，经过前向分解和后向综合。第一级前向分解与多级维纳滤波器一样；第一级之后的分解是对观测数据X_i-1(k)进行分解得到d_i(k)，d_i(k)进一步和第i级匹配滤波器h_i相乘，做减法运算，得到观测数据X_i(k)。后向综合过程与图1中一样，由一组递推的标量综合得到维纳滤波器的输出误差信号ε_i，最后一级的观测数据与权系数w_r相乘并与参考信号做差值得到误差信号ε_r-1，依次迭代，最终得到系统输出误差值ε₀在CSA-MWF结构中取r级截断，即使得X_r-1＝d_r＝ε_r。其中，i＝1,2,…,r≤N。

本发明是一种基于变对角加载量的空时抗干扰方法，流程如图3所示，包括以下几个步骤：

步骤一：通过天线对GPS卫星信号进行采样，得到采样数据，输入至相关相减多级维纳滤波器；

在空时处理中，设有M个阵元天线，每个阵元设有L个延时抽头系数，则共产生M×L-1个自由度。GPS卫星信号在奈奎斯特采样速率下进行，在时刻k各个天线的采样构成采样数据矢量，为：

$X\left(k\right)=\mathrm{AS}\left(k\right)+N\left(k\right)=a\left({\mathrm{\θ}}_0\right)s_0\left(k\right)+\underset{j=1}{\overset{L}{\mathrm{\Σ}}}a\left({\mathrm{\θ}}_j\right)s_j\left(k\right)+N\left(k\right)---\left(1\right)$

其中：X(k)为时刻k各个天线数据矢量，A＝[a(θ₀),a(θ₁),…,a(θ_L)]为阵列流形矩阵，a(θ₀)为信号导向矢量，为第j个干扰信源的导向矢量，且λ为信号波长，d为阵元间距，θ_j为第j个干扰入射角度；S(k)为信号向量，s₀(k)为期望信号向量，s_j(k)为第j个干扰信号向量，j＝1,2,…,L；N(k)为噪声向量。

线性约束中期望信号的方向和对应的导向矢量a(θ₀)认为是已知的。这里采用一维线阵，只考虑方位角。

采样数据X(k)输入到相关相减多级维纳滤波器中。

步骤二：对采样数据进行相关相减多级维纳滤波器的前向分解，得到每级的数据信号和参考信号；

现有CSA-MWF结构中阻塞矩阵为N维方阵，各级数据X_i(k)也是N维方阵，因此没有达到维数降低的目的。本发明对CSA-MWF进行优化，使得阻塞矩阵B_i为(N-i-1)×(N-i)维方阵，这样既有CSA-MWF结构的优点，又会使得各级数据X_i逐次降低，计算量减少，存储空间缩减。

具体为：

（1）设定第i级阻塞矩阵B_i；

在多级维纳滤波器中，经过信号的分解和迭代运算后，信号维数下降，运算简化。但在CSA-MWF结构中，第i级阻塞矩阵B_i为其中h_i为归一化的期望信号的导向矢量。阻塞矩阵为N维方阵，因此信号的维数并不会随着信号的分解和迭代而降低。这里取阻塞矩阵B_i维数为(N-i-1)×(N-i)，进行简化。具体处理为：取阻塞矩阵B₀前N-1行，(N-1)×N的矩阵代替原来的N维方阵；B_i用(N-i-1)×(N-i)为矩阵代替原来(N-i)×(N-i)维方阵。通过计算，这种替代保留了数据X(k)的全部信息，没有造成损失，因此是合理的。

（2）进行前向分解，得到每级的数据信号和参考信号；

前向分解过程是利用已知的期望信号信息把阵列接收到的信号分为两个支路，即上支路和下支路。其中，上支路信号X(k)经过变换后得到参考信号d₀(k)，d₀(k)包含期望信号和干扰，下支路通过阻塞矩阵B₀阻塞掉期望信号，则X₀(k)只含有干扰。如图1所示。

多级维纳滤波器中，参考信号数据信号X₀(k)＝B₀X(k)。式中为归一化的期望信号的导向矢量；B₀为满秩阻塞矩阵；||||表示向量求模。

利用相关相减算法实现前向分解如图2。CSA-MWF的阻塞矩阵为根据递推关系有：

$X_i\left(k\right)=B_i{X}_{i-1}\left(k\right)=(I-h_i{h}_i^H)X_{i-1}\left(k\right)=X_{i-1}\left(k\right)-h_i{d}_i\left(k\right)---\left(2\right)$

由上式可以看出相关相减多级维纳滤波器只需要进行第一级阻塞矩阵的计算就可以进行数据的前向分解。

根据对阻塞矩阵的改进，得到

$X_i\left(k\right)=B_i{X}_{i-1}\left(k\right)=X_{i-1}^{N-i-1}\left(k\right)-h_i^{N-i-1}{d}_i\left(k\right)---\left(3\right)$

其中，上标N-i-1表示取矩阵的前N-i-1行。

综合以上分析，前向分解为：初始条件为， X(k)；迭代运算为，X_i(k)＝X_i-1(k)-h_id_i(k)；在最后一级截断时，ε_r(k)＝d_r(k)。

其中，是前一级参考信号和观测数据的互相关函数归一化矢量，为参考信号与数据信号的自相关矩阵，i＝1,2…,r≤N-1。

步骤三：获取采样信号的协方差矩阵，采用采样信号的协方差矩阵代替真实协方差矩阵；

对角加载是为了对协方差矩阵估计进行修正，实现误差加载算法，即在降维后的协方差矩阵再加上一个对角阵。R为协方差矩阵真实量，工程中，该值是无法已知。在实际中R通过采样信号的协方差矩阵的对角元素的均值进行估计。因为协方差矩阵的对角元素是每一个阵元信号的自相关值，且每一个元素都存在误差，误差矩阵的均值为0。

（1）真实协方差矩阵的估计值采用采样信号协方差矩阵对角元素的平均值进行估计：

$\stackrel{\‾}{R}(i,i)=\mathrm{trace}\left(R\right)/\mathrm{NL}$

其中：为采样信号协方差矩阵对角元素的平均值，trace(R)为R的对角元素。

（2）判断是否对采样信号协方差矩阵进行对角加载；

如果在常规高速采样条件下，数据样本足够多，即采样频率足够大时，采样协方差矩阵能够近似为真实协方差矩阵，不需要进行对角加载，但是在小快拍下，样本数据变少，使得数据协方差矩阵的大特征值和小特征值差距很大，这就需要进行对矩阵的处理。判断方法为：

当采样频率低于信号最高频率的2倍时，信号失真，需要进行小快拍处理（例如，如信号的最高频率为1千赫兹，当采样频率低于2千赫兹时就需要进行小块拍处理），进入步骤四，获取加载量，对采样信号协方差矩阵进行对角加载。否则，直接进入步骤五，采样信号协方差矩阵的加载量设为0。

本发明基于小快拍条件下，在对卫星信号进行采样时，即计算得到采样信号协方差矩阵。

步骤四：通过矩阵求逆方法，获取对角加载的范围，进而利用干扰信号和噪声信号特征值的关系自适应的选择需要加载的对角量；

变对角加载量是根据特征值的大小选择加载量可以进一步减少噪声的扩散度，降低加载量对大特征值的影响。

步骤4.1：确定加载量范围

对角加载后矩阵R_dl：

R_dl＝R+εB+γI(4)

其中，R为协方差矩阵真实量；ε是采样协方差矩阵误差常量，B是均值为0、方差为1的随机矩阵；γ是对角加载方法的加载量，I为单位对角阵。固定加载量是通过经验值进行加载，典型值是背景噪声的5～10倍即取这种选择既保证了干扰信号的大特征值不受影响，又能够降低噪声信号对小特征值的发散度，提高了滤波器性能。但是固定加载值是根据经验值进行确定，在实际工程中会出现一些噪声和复杂干扰的情况，固定加载量并不能完全有效的改善干扰影响。针对这种问题，本发明提出可变对角加载方法。

通过矩阵求逆定理，对(4)式进行求逆可得到对角加载协方差矩阵的近似表达式：

$R_{\mathrm{dl}}^{-1}=(R+\mathrm{\γI})\{I-\frac{\mathrm{\ϵ}}{\mathrm{\γ}+{\mathrm{\σ}}_n^2}B[I-V{(V^{H}V+(\mathrm{\γ}+{\mathrm{\σ}}_n^2){\mathrm{\Λ}}^{-1})}^{-1}{V}^H\left]\right\}---\left(5\right)$

其中，V和Λ是矩阵分解过程中的特征向量，V为信号和干扰子空间特征向量，Λ为在噪声子空间特征向量。对角加载值应小于协方差矩阵的对角元素，即γ＜R(i,i),i＝i,2,…,NL。由（5）式可以看出，滤波器性能的降低由大括号内的多项式引起，如果多项式第二项为0，则该表达为最优波束形成。因此得到对角加载的范围ε≤γ＜R(i,i),i＝1,2,…,NL。

步骤4.2：对加载后协方差矩阵对角化

CSA-MWF的阻塞矩阵选择所以显然有[h₁,h₂,…,h_r]是标准正交向量组。从而可以用矩阵T_D对（4）式对角加载后的协方差矩阵R_dl进行三角化

$T_D^H{R}_{\mathrm{dl}}{T}_D=T_D^H(R_{X_0}+\mathrm{\γ}{I}_{(M-1)L+1})T_D=R_d+{\mathrm{\γI}}_r$

式中，T_D为变换矩阵；为第i级参考信号的自相关量，δ_i为第i+1级输出误差值与第i级参考信号相关量，i＝1,2,…,r；^*为求矩阵共轭。进一步用变对角加载函数表示：

$T_D^H{R}_{\mathrm{dl}}{T}_D=T_D^H(R_{X_0}+T_D{\mathrm{BT}}_D^H)T_D=R_d+B$

其中，矩阵f(.)为变加载量函数。

步骤4.3：计算加载函数

加载函数：

$f\left({\mathrm{\σ}}_{d_i}^2\right)=\left\{\begin{array}{cc}0& \mathrm{else}\\ \mathrm{\γ}-({\mathrm{\σ}}_{d_i}^2-{\mathrm{\σ}}_{d_r}^2)& \mathrm{\γ}-({\mathrm{\σ}}_{d_i}^2-{\mathrm{\σ}}_{d_r}^2)\end{array}\right.---\left(8\right)$

对角加载后多级维纳滤波器输出结果ε_i方差表达为

$E\left[{\left|{\mathrm{\ϵ}}_i\left(k\right)\right|}^2\right]=\frac{1}{K}{\mathrm{\ϵ}}_i{{\mathrm{\ϵ}}_i}^H=\frac{1}{K}{\mathrm{\ϵ}}_i{{\mathrm{\ϵ}}_i}^H+{\mathrm{\γ}}_i=E\left[{\left|{\mathrm{\ϵ}}_i\left(k\right)\right|}^2\right]+{\mathrm{\γ}}_i---\left(9\right)$

式中，ε_i为维纳滤波器第i级输出误差信号，K为信号采样快拍数，表达式 ${\mathrm{\γ}}_i=f\left({\mathrm{\σ}}_{d_i}^2\right)+f\left({\mathrm{\σ}}_{d_{i+1}}^2\right){\left|w_{i+1}\right|}^2+...+f\left({\mathrm{\σ}}_{d_r}^2\right)\underset{m=i+1}{\overset{r}{\mathrm{\Π}}}{\left|w_m\right|}^2.$

由此可以得到递推关系：

$\begin{array}{c}{\mathrm{\γ}}_r=f\left({\mathrm{\σ}}_{d_r}^2\right)\\ {\mathrm{\γ}}_i=f\left({\mathrm{\σ}}_{d_i}^2\right)+{\mathrm{\γ}}_{i+1}{\left|w_{i+1}\right|}^2\end{array}---\left(10\right)$

其中，w_i+1为每一级滤波器的自适应权值；

步骤4.4：对加载量检测

由9式、10式做迭代运算，得到CSA-MWF加载矩阵之后，依据步骤4.1检验是否在加载范围内，如果在加载范围内，则数据有效；如果数据较ε小，则数据直接取ε；如果数据较R(i,i)大，则数据直接取R(i,i)。

此过程中得到随干扰空间和噪声空间特征值变化的对角元素的加载量，并加载到步骤三中的采样协方差矩阵上，得到加载后的新协方差矩阵。

步骤五：进行相关相减多级维纳滤波器的后向综合，求出滤波器最优权矢量，对当前信号进行自适应滤波；

计算权CSA-MWF的矢量：

$w_i=R_{{\mathrm{\ϵ}}_i}^{-1}{r}_{{\mathrm{\ϵ}}_i{d}_{i-1}}={\mathrm{\ξ}}_i^{-1}{\mathrm{\δ}}_i---\left(11\right)$

$r_{{\mathrm{\ϵ}}_{i+1}{d}_i}=E\left[{\mathrm{\ϵ}}_{i+1}\left(k\right)d_i^{*}\left(k\right)\right]=h_{i+1}^H{r}_{X_i{d}_i}=\sqrt{r_{X_i{d}_i}^H{r}_{X_i{d}_i}}={\mathrm{\δ}}_{i+1}---\left(12\right)$

其中，w_i为滤波器的第i级权系数，所以，

$\begin{array}{c}{\mathrm{\ξ}}_i=E\left[{\left|{\mathrm{\ϵ}}_i\left(k\right)\right|}^2\right]=E\left[\right|d_i\left(k\right)-w_{i+1}^{*}{\mathrm{\ϵ}}_{i+1}\left(k\right){|}^2]\\ ={\mathrm{\σ}}_{d_i}^2-w_{i+1}^{*}{\mathrm{\δ}}_{i+1}={\mathrm{\σ}}_{d_i}^2-{\mathrm{\ξ}}_{i+1}^{-1}{\left|{\mathrm{\δ}}_{i+1}\right|}^2\end{array}----\left(13\right)$

其中，ξ_i为第i级输出误差的自相关矩阵，系统每级输出值这里做r级切断，即ε_r(k)＝d_r(k)。r＝N-1时，为满秩处理。滤波器结构的输出ε₀可以由递推公式得到

${\mathrm{\ϵ}}_0=d_0\left(k\right)-w_1^{*}{\mathrm{\ϵ}}_1\left(k\right)=d_0\left(k\right)-w_1^{*}(d_1\left(k\right)-w_2^{*}{\mathrm{\ϵ}}_2\left(k\right))$

=...

$=d_0\left(k\right)-w_1^{*}\{d_1\left(k\right)-w_2^{*}[d_2\left(k\right)-...w_{r-1}^{*}(d_{r-1}\left(k\right)-w_r^{*}{\mathrm{\ϵ}}_r\left(k\right))\}---\left(14\right)$

$=d_0\left(k\right)-W_d^{H}d\left(k\right)$

其中，W_d是关于滤波器权系数的矩阵，表达式为d(k)是关于参考信号的矩阵，表达式为d(k)＝[d₁(k),d₂(k),…,d_N(k)]^T,第i级的参考信号为

$d_i\left(k\right)=h_i^H{X}_{i-1}\left(k\right)=h_i^H{B}_{i-1}{X}_{i-2}\left(k\right)=...=h_i^H\left(\underset{j=i-1}{\overset{1}{\mathrm{\Π}}}{B}_j\right)X_0\left(k\right).$

所以，

$d\left(k\right)={[h_1^H,h_2^H{B}_1,...,h_{n-1}^H\underset{j=n-2}{\overset{1}{\mathrm{\Π}}}{B}_j,\underset{j=n-1}{\overset{1}{\mathrm{\Π}}}{B}_j]}^T{X}_0\left(k\right)=L{X}_0\left(k\right)---\left(15\right)$

引入步骤四中求得的对角加量后，后向递推过程中wi的关系式改变：

$w_i=R_{{\mathrm{\ϵ}}_i}^{-1}{r}_{{\mathrm{\ϵ}}_i{d}_{i-1}}={\mathrm{\ξ}}_i^{-1}{\mathrm{\δ}}_i$

$=E\left[d_{i-1}^{*}\left(k\right){\mathrm{\ϵ}}_i\left(k\right)\right]/E\left[{\left|{\mathrm{\ϵ}}_i\left(k\right)\right|}^2\right]$

$=E\left[d_{i-1}^{*}\left(k\right)(d_{i-1}\left(k\right)-w_i^{*}{\mathrm{\ϵ}}_i\left(k\right))\right]/\left(E\right[{\left|{\mathrm{\ϵ}}_{i+1}\left(k\right)\right|}^2]+{\mathrm{\γ}}_i)---\left(16\right)$

$=E\left[d_{i+1}^{*}\left(k\right)d_i\left(k\right)\right]/\left(E\right[{\left|{\mathrm{\ϵ}}_i\left(k\right)\right|}^2]+{\mathrm{\γ}}_i)$

＝δ_i/(E[|ε_i(k)|²]+γ_i)

由此，可以得到CSA-MWF结构的权矢量

$W_{\mathrm{MWF}}=h_0-B_0^H{L}^H{W}_d$

其中，符号L为 $L={[h_1^H,h_2^H,B_1,...,h_{n-1}^H\underset{j=n-2}{\overset{1}{\mathrm{\Π}}}{B}_j,\underset{j=n-1}{\overset{1}{\mathrm{\Π}}}{B}_j]}^T,$ Π表示连乘符号。

权矢量W_MWF体现滤波器的滤波参数，当系统输出的误差值ε₀符合某种条件时，该滤波器的性能最优。本发明是在线性约束最小方差（LCMV）准则下，使ε₀最小即为最优滤波，即要求天线阵列在任何情况下要确保信号方向的增益为1条件下，使ε₀最小。

空时抗干扰是在空域和时域对采样信号进行滤波。在采样点较小或者很小时，采样得到的信息不足以使得滤波器结构达到最优，本发明引入对角加载方法弥补小快拍下信息不足带来的影响。

Claims (1)

1.一种基于变对角加载量的空时抗干扰方法，其特征在于，包括以下几个步骤：

步骤一：通过天线对GPS卫星信号进行采样，得到采样数据，输入至相关相减多级维纳滤波器；

设有M个阵元天线，每个阵元设有L个延时抽头系数，则共产生M×L-1个自由度；GPS卫星信号在奈奎斯特采样速率下进行，在时刻k各个天线的采样构成采样数据矢量，为：

其中：X(k)为时刻k各个天线数据矢量，A＝[a(θ₀),a(θ₁),…,a(θ_L)]为阵列流形矩阵，a(θ₀)为信号导向矢量，为第j个干扰信源的导向矢量，且λ为信号波长，d为阵元间距，θ_j为第j个干扰入射角度；S(k)为信号向量，s₀(k)为期望信号向量，s_j(k)为第j个干扰信号向量，j＝1,2,…,L；N(k)为噪声向量；

设定线性约束中期望信号向量的方向和对应的导向矢量a(θ₀)为已知，采用一维线阵，只考虑方位角；

采样数据矢量X(k)输入到相关相减多级维纳滤波器中；

步骤二：对采样数据进行相关相减多级维纳滤波器的前向分解，得到每级的数据信号和参考信号；

具体为：

(1)设定第i级阻塞矩阵B_i；

第i级阻塞矩阵B_i为(N-i)×(N-i)维方阵，取其(N-i-1)×(N-i)维，设为新的第i级阻塞矩阵B_i；

(2)进行前向分解，得到每级的数据信号和参考信号；

前向分解是利用已知的期望信号信息把阵列接收到的信号分为两个支路，即上支路信号和下支路信号；其中，上支路信号为采样数据矢量X(k)，X(k)经过变换后得到参考信号d₀(k)，下支路信号通过阻塞矩阵B₀阻塞掉期望信号得到数据信号X₀(k)；

参考信号数据信号X₀(k)＝B₀X(k)；式中h₀为归一化的期望信号的导向矢量；B₀为满秩阻塞矩阵；||||表示向量求模；

利用相关相减算法实现前向分解，相关相减多级维纳滤波器的阻塞矩阵为根据递推关系有：

根据对阻塞矩阵的改进，得到

其中，上标N-i-1表示取矩阵的前N-i-1行；

前向分解为：初始条件为， X₀(k)＝B₀X(k)；迭代运算为，X_i(k)＝X_i-1(k)-h_id_i(k)；在最后一级截断时，ε_r(k)＝d_r(k)；

其中，是前一级参考信号和观测数据的互相关函数归一化矢量，i＝1,2…,r≤N-1；

步骤三：获取采样信号的协方差矩阵，采用采样信号的协方差矩阵代替真实协方差矩阵；

具体为：

(1)真实协方差矩阵的估计值采用采样信号协方差矩阵对角元素的平均值进行估计：

其中：为采样信号协方差矩阵对角元素的平均值，trace(R)为R的对角元素；

(2)判断是否对采样信号协方差矩阵进行对角加载；

当采样频率低于信号最高频率的2倍时，进入步骤四，获取加载量，对采样信号协方差矩阵进行对角加载；否则，直接进入步骤五，采样信号协方差矩阵的加载量设为0；

步骤四：通过矩阵求逆方法，获取对角加载的范围，进而利用干扰信号和噪声信号特征值的关系自适应的选择加载量；

步骤4.1：确定加载量范围；

加载后的采样信号协方差矩阵R_dl：

R_dl＝R+εB₁'+γI(4)

其中，R为采样信号协方差矩阵；ε是采样信号协方差矩阵误差常量，B₁'是均值为0、方差为1的随机矩阵；γ是对角加载方法的加载量，I为单位对角阵；

通过矩阵求逆定理，对(4)式进行求逆，得到加载后的采样信号协方差矩阵的近似表达式：

其中，V和Λ是矩阵分解过程中的特征向量，V为信号和干扰子空间特征向量，Λ为在噪声子空间特征向量，加载量应小于采样信号协方差矩阵的对角元素，即γ＜R(i,i),i＝ 1,2,…NL，得到对角加载的范围ε≤γ＜R(i,i),i＝1,2,…,NL；

步骤4.2：对加载后协方差矩阵对角化；

采用矩阵T_D对(4)式对角加载后的采样信号协方差矩阵R_dl进行三角化

式中，T_D为变换矩阵；为第i级参考信号的自相关量，δ_i为第i+1级输出误差值与第i级参考信号相关量，i＝1,2,…,r；^*为求矩阵共轭；进一步用加载函数表示：

其中，矩阵f(.)为加载函数；

步骤4.3：计算加载函数；

加载函数：

对角加载后多级维纳滤波器第i级输出误差信号方差表达为

式中，为对角加载后多级维纳滤波器第i级输出误差信号，K为信号采样快拍数，

得到递推关系：

其中，w_i+1为滤波器的第i+1级权系数；

步骤4.4：对加载量检测

由式(9)、(10)做迭代运算，得到相关相减多级维纳滤波器加载矩阵之后，依据步骤4.1检验是否在加载范围内，如果在加载范围内，则数据有效；如果数据较ε小，则数据直接取ε；如果数据较R(i,i)大，则数据直接取R(i,i)；

将得到随干扰信号和噪声信号特征值变化的对角元素的加载量，加载到步骤三中的采样信号协方差矩阵上，得到加载后的新协方差矩阵；

步骤五：进行相关相减多级维纳滤波器的后向综合，求出滤波器最优权矢量，对当前信号进行自适应滤波；

计算相关相减多级维纳滤波器结构的权矢量：

其中，w_i为滤波器的第i级权系数，所以，

其中，ξ_i为第i级输出误差值的自相关矩阵，系统每级输出误差值做r级切断，即ε_r(k)＝d_r(k)；r＝N-1时，为满秩处理，滤波器结构的输出ε₀由递推公式得到：

其中，W_d是关于滤波器权系数的矩阵，表达式为是关于参考信号的矩阵，表达式为d(k)＝[d₁(k),d₂(k),…,d_N(k)]^T,第i级的参考信号为

所以，

引入加载后的采样信号协方差矩阵，后向递推过程中w_i的关系式改变：

由此，得到相关相减多级维纳滤波器结构的权矢量：

其中，符号L为∏表示连乘符号。