CN104155632A - 一种基于局部相关性的改进子空间海杂波抑制方法 - Google Patents

Abstract

本发明公开一种基于局部相关性的改进子空间海杂波抑制方法。海杂波的抑制能够提高目标信杂比，有助于提高高频地波雷达(high frequency surface wave radar,HFSWR)对舰船目标的检测性能。原始子空间法抑制海杂波时对所有距离单元构造的协方差矩阵是相同的，造成对杂波空间估计不准从而抑制效果不佳，距离多普勒谱上反映为不同距离海杂波抑制效果类似。本发明通过计算待检测距离单元与相邻距离单元之间的相关系数，自适应的确定组成协方差矩阵的参考单元，对构成的协方差矩阵进行奇异值分解(singular value decomposition,SVD)分解得到海杂波所在的子空间，并最终将其抑制。该方法具有同时抑制一阶和高阶海杂波的能力，提高了目标信杂比，对所有距离单元的抑制效果都较理想。

Description

一种基于局部相关性的改进子空间海杂波抑制方法

技术领域：

本发明涉及一种高频地波雷达海杂波抑制方法，特别涉及一种基于局部相关性的改进子空间海杂波抑制方法。

背景技术：

高频地波雷达工作在3～30MHz频段，是超视距雷达的一种，其利用垂直极化的电磁波沿海表面绕射传播损耗小的特点，可实现对海洋动力学参数、海上慢速目标及低空飞机的监视与定位，近年来广泛应用于军事与民用领域。

高频地波雷达回波情况复杂，除了包括射频干扰、瞬态干扰和电离层干扰等各类干扰外，还包括很强的海杂波。海杂波能量主要集中在两个由雷达工作频率决定的一阶海杂波，即Bragg峰处，实际情况时，一阶Bragg峰会有不同程度的展宽。对于舰船目标检测来说，由于舰船目标的多普勒频率和强大的海杂波相近，海杂波才是限制检测性能的主要因素。因此相对于飞机目标，舰船目标检测也更为困难。所以，在进行舰船目标检测时，一个重要的思想就是要先抑制海杂波，提高信杂比，以期对后面的检测有所帮助。

现有较常用的海杂波抑制方法主要有以下两个思路：一类是时域对消类方法，如循环对消方法和基于幅度相位估计(amplitude and phase estimation,APES)的对消方法。这类方法通过估计强杂波对应谐波分量的频率、幅度和相位3个参数，有选择性的对消掉强杂波分量，使目标信号从杂波中显露出来，这类方法的性能取决于上述3个参数的估计精度。(如文献1：Root B.HF radar ship detection through clutter cancellation[C]//RadarConference,1998.RADARCON98.Proceedings of the1998IEEE.IEEE,1998:281-286.)(如文献2：郭欣,倪晋麟,刘国岁.短相干积累条件下天波超视距雷达的舰船检测[J].电子与信息学报,2004,26(4):613-618)。另一类是基于矩阵分解的海杂波方法，如基于子空间的方法、基于奇异值分解(singular-value decomposition,SVD)的方法和通过多普勒频率估计杂波子空间(estimate subspace via Doppler,ESVID)的方法，这类方法根据杂波在子空间的聚集特性来实现海杂波抑制。其中，基于子空间方法通过参考单元数据来估计杂波子空间，(如文献3：邢孟道,保铮,强勇.天波超视距雷达瞬态干扰抑制[J].电子学报,2002,30(6):823-826.)；SVD方法通过待检测单元来估计杂波子空间，(如文献4：Poon M W Y,Khan R H,Le-Ngoc S.A singular value decomposition(SVD)basedmethod for suppressing ocean clutter in high frequency radar[J].Signal Processing,IEEETransactions on,1993,41(3):1421-1425.)；ESVID方法根据杂波子空间与信号矢量张成子空间的等价性，通过多普勒频率估计来估计杂波子空间，(如文献5：Zhiguo Z,Chen J,Zheng B.A method to estimate subspace via Doppler for ocean clutter suppression inskywave radars[C]//Radar(Radar),2011IEEE CIE International Conference on.IEEE,2011,1:145-148.)。

除了上述较常用的两大类方法外，一些新颖的方法也被提出。如中国专利公开CN102928823B基于二次雷达应答信号相关特征的海杂波抑制方法中提出将受到海杂波影响的探测目标计算其真实应答信号和受海杂波效应的反射信号之间的相对时间差，根据不同相对时间差的码位值演变方式，生成与真实目标相对应的目标航管代码并保存；将真实目标航管代码与历史探测点迹数据进行代码相关，相关成功后对探测目标航管代码进行纠码并实现海杂波的抑制。如中国专利公开CN103645467A海杂波抑制以及海杂波背景中目标检测的方法和系统中提出对海面搜索雷达实测的一维距离像的历程图进行二维傅里叶变换后得到估计的海杂波一维距离像的历程图；将实测的一维距离像的历程图的图像数据减去估计的海杂波一维距离像的历程图的图像数据，得到抑制海杂波的一维距离像的历程图的图像数据。如中国专利公开CN102779336A基于混沌的光电图像海杂波抑制方法中提出利用光电图像海杂波所具有的混沌特性，采用RBF神经网络方法预测与抑制海杂波。如中国专利公开CN102169175B基于FRFT域非相参积累的海杂波抑制和动目标检测方法中通过对FRFT域延迟数据的迭代平滑，分别采用FRFT域均值积累方法和FRFT域单极点反馈积累方法，完成非相参积累，将FRFT域信号与积累后的数据对消，抑制海杂波。

时域对消类海杂波抑制方法的效果很难保证，当脉冲相干积累时间较短时则频率分辨力低，此时对Bragg峰的频率和幅度进行估计时精度下降，无法得到较好的初相估计，造成循环对消次数增加，而随着对消的不断进行，剩余杂波将会在多普勒中扩散开来，影响目标的检测。

子空间的海杂波抑制方法中，对于一个进行完距离向脉冲压缩和方位向波束形成后的三维数据块，取整个数据块进行协方差估计，即认为整个数据块的海杂波是高度相关的，用的是同一个标准进行海杂波抑制，这显然是不合理的。海杂波确实在一定的范围内高度相关，但这个范围并没有大到整个雷达探测的范围。对于整体用一个标准来进行海杂波抑制，有可能出现近处海杂波抑制不干净或远处海杂波抑制过多的情况。

那么，对于不同距离不同方位的待抑制回波数据，应该设计不同的参考单元个数来构造与其对应的协方差矩阵，从而进行海杂波的抑制。参考单元的个数，以及协方差矩阵的构造直接影响到海杂波抑制效果。参考单元个数取的过大，海杂波空间相关性减弱，协方差估计不准，在很大一片区域都用同一个标准进行抑制不符合海杂波在空间分布的特点；参考单元个数取的过小，则样本数太小也无法进行准确估计。

发明内容：

本发明为了提高高频地波雷达舰船目标检测性能，提出了基于局部相关性的改进子空间海杂波抑制方法，对于每个特定的方位和距离单元的待抑制数据，通过计算周围单元与待抑制单元的相关系数，自适应的调整块的大小以期获得更精确的协方差矩阵估计，继而进行后续的海杂波抑制。

本发明采用如下技术方案：一种基于局部相关性的改进子空间海杂波抑制方法，其包括如下步骤：

(a)、设输入数据为距离向脉冲压缩，方位向数字波束形成后的距离-慢时间-方位三维数据块；

(b)、每次抑制的对象为一固定方位固定距离门的一维慢时间数据，

设当前待检测的某确定方位和距离的距离单元的K个回波慢时间数据为矢量S_c＝[s_c(0),s_c(1),...,s_c(K-1)]^T，第i个距离单元的K个回波慢时间数据为S_i＝[s_i(0),s_i(1),...,s_i(K-1)]^T，其中T为转置；

(c)、遍历出与当前待检测距离单元相邻的且方位与待检测距离单元相同的距离单元慢时间数据S_i，并计算其与待检测距离单元S_c的相关系数

${\mathrm{\ρ}}_i=\frac{{S_c}^H{S}_i}{{\left|\right|S_c\left|\right|}_2\·{\left|\right|S_i\left|\right|}_2}$

其中H代表共轭转置，||||₂代表二范数，其中0≤|ρ_i|≤1，1表示相关性最大，0表示相关性最小，设计门限a找到ρ≥a的P个对应S_c的参考单元数据S_i(i＝1,2,...P)，和对应的相关系数ρ_i(i＝1,2,...P)；

(d)、按距离从近到远将所选P个参考单元数据排列，得到X＝[S₁,S₂,...S_p]，按照对应的顺序，构造W＝diag(ρ₁,ρ₂,...ρ_P)，则S_c的估计协方差矩阵为

(e)、将R进行SVD奇异值分解

R＝VΣV^H

Σ＝diag(λ₁,λ₂,...,λ_K)，得到K个奇异值(λ₁≥λ₂≥...≥λ_K)，V＝[v₁,v₂,...v_K]，v_i为对应λ_i的特征向量，且V为酉矩阵，所述SVD分解为矩阵的一种正交分解，矩阵的主要分量与大奇异值相对应，定义归一化奇异值为：

$\stackrel{\‾}{{\mathrm{\λ}}_i}=\frac{{\mathrm{\λ}}_i}{{\mathrm{\λ}}_1}$

则，选择一个小正数作为阈值σ，则比σ大的前r个特征值对应的特征向量所张成的空间为信号子空间即海杂波所在空间V_c＝[v₁,v₂,...v_r]，其余特征向量所张成的空间为噪声子空间；

(f)、将待抑制的S_c投影到信号子空间V_c＝[v₁,v₂,...v_r]上，得到滤波后的输出为

S_{c_sup}＝(I-V_cV_c ^H)S_c

其中I为K×K维单位矩阵；

(g)、重复(b)～(f)的步骤，检测当前距离-慢时间-方位三维数据块中的待检测方位上的所有距离单元数据，接着选择其它方位，直至得出最后一个方位上的最后一个距离单元的数据，则当前距离-慢时间-方位三维数据块海杂波抑制完成。

本发明具有如下有益效果：

(1).算法可以实现高频地波雷达海杂波抑制，提高了信杂比，有助于后续舰船目标检测，提高了系统的鲁棒性；

(2).由于针对每个方位和距离单元自适应的调整协方差估计的参数，更为符合海杂波在空间上的分布特征，相比原来整个数据块进行抑制时效果更好。

附图说明：

图1高频地波雷达舰船目标检测流程图。

图2待进行海杂波抑制的时域三维数据块。

图3待抑制单元与参考单元位置示意图。

图4基于局部相关性的改进子空间海杂波抑制方法流程图。

图5当前待抑制距离单元慢时间数据与各距离单元慢时间数据相关系数。

图6奇异值分解后的归一化奇异值分布情况。

图7未抑制海杂波时的时域谱。

图8抑制完海杂波后的时域谱。

图9海杂波抑制前后doppler谱对比图。

图10改进方法海杂波抑制前后距离doppler谱图对比。

图11原始子空间方法海杂波抑制前后距离doppler谱图对比。

具体实施方式：

下面结合附图对本发明的技术方案进行详细说明：

图1中可以揭示出海杂波抑制模块在高频地波雷达舰船目标检测信号处理流程中的位置，海杂波抑制处理时在距离向脉冲压缩，方位向数字波束形成之后的操作。将一帧数据(也称为一个驻留时间)定义为快时间-慢时间-通道的三维数据块。其中快时间采样个数对应距离单元的个数，慢时间采样个数对应雷达相干积累的脉冲个数，通道数为雷达接收天线的个数。本发明提出的海杂波抑制方法是时域的海杂波抑制方法，即在快时间维脉冲压缩得到距离信息，通道维数字波束形成得到方位信息之后进行海杂波的抑制。待进行海杂波抑制的时域三维数据块如图2所示。抑制完海杂波后，再进行慢时间维多普勒处理得到速度信息，之后进行恒虚警检测、点迹融合、航迹关联等目标检测算法，最终输出目标。

本发明基于局部相关性的改进子空间海杂波抑制方法，其包括如下步骤：

(a)、设输入数据为距离向脉冲压缩，方位向数字波束形成(digital beam forming,DBF)后的距离-慢时间-方位三维数据块，如图2所示。

(b)、每次抑制的对象为一固定方位固定距离门的一维慢时间数据。

设当前待检测的某确定方位和距离的距离单元的K个回波慢时间数据为矢量S_c＝[s_c(0),s_c(1),...,s_c(K-1)]^T，第i个距离单元的K个回波慢时间数据为S_i＝[s_i(0),s_i(1),...,s_i(K-1)]^T，其中T为转置。

(c)、遍历出与当前待检测距离单元相邻的且方位与待检测距离单元相同的距离单元慢时间数据S_i，并计算其与待检测距离单元S_c的相关系数

${\mathrm{\ρ}}_i=\frac{{S_c}^H{S}_i}{{\left|\right|S_c\left|\right|}_2\·{\left|\right|S_i\left|\right|}_2}$

其中H代表共轭转置，||||₂代表二范数。0≤|ρ_i|≤1，1表示相关性最大，0表示相关性最小，设计门限a找到ρ≥a的P个对应S_c的参考单元数据S_i(i＝1,2,...P)，和对应的相关系数ρ_i(i＝1,2,...P)。待抑制单元与参考单元的结构示意图如图3所示。

(d)、按距离从近到远将所选P个参考单元数据排列，得到X＝[S₁,S₂,...S_p]，按照对应的顺序，构造W＝diag(ρ₁,ρ₂,...ρ_P)，则S_c的估计协方差矩阵为

(e)、将R进行SVD(奇异值)分解

R＝VΣV^H

Σ＝diag(λ₁,λ₂,...,λ_K)，得到K个奇异值(λ₁≥λ₂≥...≥λ_K)。V＝[v₁,v₂,...v_K]，v_i为对应λ_i的特征向量，且V为酉矩阵。SVD分解实际上是矩阵的一种正交分解，矩阵的主要分量与大奇异值相对应。定义归一化奇异值为：

$\stackrel{\‾}{{\mathrm{\λ}}_i}=\frac{{\mathrm{\λ}}_i}{{\mathrm{\λ}}_1}$

显然，选择一个小正数作为阈值σ(例如0.05)，则认为比σ大的前r个特征值对应的特征向量所张成的空间为信号子空间(海杂波所在空间V_c＝[v₁,v₂,...v_r]，其余特征向量所张成的空间为噪声子空间(目标和噪声所在空间)。

(f)、将待抑制的S_c投影到信号子空间V_c＝[v₁,v₂,...v_r]上，得到滤波后的输出为

S_{c_sup}＝(I-V_cV_c ^H)S_c

其中I为K×K维单位矩阵。

(g)、重复(b)～(f)的步骤，检测当前距离-慢时间-方位三维数据块中的待检测方位上的所有距离单元数据，接着选择其它方位，直至得出最后一个方位上的最后一个距离单元的数据，则当前距离-慢时间-方位三维数据块海杂波抑制完成。

整个海杂波抑制流程如图4所示。

取某高频地波雷达实测数据进行算法检验。该高频地波雷达共有80个距离门，距离分辨率为5km，则距离向探测范围为0～400km；1个发射天线12根接收天线，13个方位单元，方位分辨率为10°，方位向探测范围为120°。依据待测海面舰船目标的径向速度范围和雷达脉冲重复频率选择相干积累脉冲数为256。所以一帧三维数据块的维数为80×256×13。

使用2013年7月12日15时17分的一帧三维数据进行基于子空间的海杂波抑制验证。根据AIS及船讯网信息已知一确定的舰船目标在(84.574km，20°)处，数据上看在第16个距离门、第152个doppler频率、第11个方位上。

取第11个方位，第16个距离门的数据作为待抑制单元进行检验，抑制前的海杂波时域谱如图7所示。计算其与周围距离单元的数据得到相关系数如图5所示，设计门限a＝0.75，得到ρ≥a的参考单元个数25个，即P＝25。进行协方差估计后，得到K个归一化奇异值分布如图6所示，取阈值σ＝0.05，则比阈值σ大的个数为r＝35，之后按前面所述步骤进行后续海杂波抑制，最终得到该抑制单元的抑制完海杂波的时域数据，如图8所示。将抑制前后的时域谱分别做快速傅里叶变换(Fast Fourier Transform，FFT)，得到doppler谱，并放在一张图中如图9所示。可以清楚的发现，一阶海杂波及其周围较高的二阶海杂波得到了很好的抑制，但是目标并没有被抑制掉，相较周围，目标被凸显出来方便了后续的检测。对第11个方位上所有的距离进行抑制，情况如图10所示。可以发现，较原始的抑制方法(效果如图11所示)，整体上海杂波得到了较好的抑制，提高了目标信杂比。

以上所述仅是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下还可以作出若干改进，这些改进也应视为本发明的保护范围。

Claims (1)

1.一种基于局部相关性的改进子空间海杂波抑制方法，其特征在于：包括如下步骤

(a)、设输入数据为距离向脉冲压缩，方位向数字波束形成后的距离-慢时间-方位三维数据块；

(b)、每次抑制的对象为一固定方位固定距离门的一维慢时间数据，

设当前待检测的某确定方位和距离的距离单元的K个回波慢时间数据为矢量S_c＝[s_c(0),s_c(1),...,s_c(K-1)]^T，第i个距离单元的K个回波慢时间数据为S_i＝[s_i(0),s_i(1),...,s_i(K-1)]^T，其中T为转置；

(c)、遍历出与当前待检测距离单元相邻的且方位与待检测距离单元相同的距离单元慢时间数据S_i，并计算其与待检测距离单元S_c的相关系数

${\mathrm{\ρ}}_i=\frac{{S_c}^H{S}_i}{{\left|\right|S_c\left|\right|}_2\·{\left|\right|S_i\left|\right|}_2}$

其中H代表共轭转置，||||₂代表二范数，其中0≤|ρ_i|≤1，1表示相关性最大，0表示相关性最小，设计门限a找到ρ≥a的P个对应S_c的参考单元数据S_i(i＝1,2,...P)，和对应的相关系数ρ_i(i＝1,2,...P)；

(d)、按距离从近到远将所选P个参考单元数据排列，得到X＝[S₁,S₂,...S_p]，按照对应的顺序，构造W＝diag(ρ₁,ρ₂,...ρ_P)，则S_c的估计协方差矩阵为

(e)、将R进行SVD奇异值分解

R＝VΣV^H

Σ＝diag(λ₁,λ₂,...,λ_K)，得到K个奇异值(λ₁≥λ₂≥...≥λ_K)，V＝[v₁,v₂,...v_K]，v_i为对应λ_i的特征向量，且V为酉矩阵，所述SVD分解为矩阵的一种正交分解，矩阵的主要分量与大奇异值相对应，定义归一化奇异值为：

$\stackrel{\‾}{{\mathrm{\λ}}_i}=\frac{{\mathrm{\λ}}_i}{{\mathrm{\λ}}_1}$

则，选择一个小正数作为阈值σ，则比σ大的前r个特征值对应的特征向量所张成的空间为信号子空间即海杂波所在空间V_c＝[v₁,v₂,...v_r]，其余特征向量所张成的空间为噪声子空间；

(f)、将待抑制的S_c投影到信号子空间V_c＝[v₁,v₂,...v_r]上，得到滤波后的输出为

S_{c_sup}＝(I-V_cV_c ^H)S_c

其中I为K×K维单位矩阵；

(g)、重复(b)～(f)的步骤，检测当前距离-慢时间-方位三维数据块中的待检测方位上的所有距离单元数据，接着选择其它方位，直至得出最后一个方位上的最后一个距离单元的数据，则当前距离-慢时间-方位三维数据块海杂波抑制完成。