CN103217670B - 一种基于pca的外辐射源微弱信号检测方法 - Google Patents

Abstract

本发明提供一种基于主量分析PCA的外辐射源微弱信号检测方法，包括以下步骤：求采样序列x(t)的自相关矩阵R=E{x(t)x^H(t)}；对自相关矩阵R进行奇异值分解SVD得到R=QΛQ^H，利用分解得到的对角阵Λ以及酉矩阵Q对采样序列x(t)进行线性变换得到接收信号从强到弱的排列y(t)；将接收信号y(t)中第3路至第M路信号分别与第1路信号y₁(t)做时频分布二维相关得到对应相关值|Ψ_i(τ,f)|，当有某一路信噪比大于等于门限值时，则将其峰值点对应的时延与多普勒频确定为一个目标的时延与多普勒频移，否则判定该路信号中不存在目标。本发明能够实现在强同频干扰以及低信噪比条件下提取微弱目标信号。

Description

一种基于PCA的外辐射源微弱信号检测方法

技术领域

本发明属于信号处理技术，具体涉及微弱信号检测技术。

背景技术

基于外辐射源的双多基地雷达***能够利用地理位置已知的民用电视信号、广播信号、通信信号、卫星信号，以及手机基站等作为照射源，通过对空中运动目标反射的回波信号的检测，实现对目标的检测与定位跟踪。由于外辐射源雷达自身不发射信号，整个雷达***具有隐蔽性强、抗干扰能力强和生存能力强等优点，因而受到各国的广泛关注。

***的原理示意如图1所示，以广播信号为例，雷达到接收信号的包括有来自于广播发射塔本身的直达波、运动目标的回波、其它同频广播发射塔发送的直达波以及其它物体反射的回波，比如各种建筑反射的静物回波。因此，雷达接收信号中的目标回波信号是极其微弱的，并且淹没于各种干扰中。

在实际应用上它面临着许多的挑战，其中，最直接的一个问题就是基于民用信号的外辐射源雷达***如何更有效地实现对微弱目标的检测与跟踪，主要体现为在外辐射源雷达的接收信号中，目标回波信号淹没在直达波和多径信号等杂波以及接收机热噪声中，增加了目标检测的难度。

传统对外辐射源目标信号的检测方法是用参考通道接收直达波信号，主天线接收多路信号，其中包括直达波信号、目标回波信号以及噪声和干扰。处理过程如图2所示，首先用主通道信号和参考通道信号采用NLMS算法对直达波进行自适应对消，再通过空域滤波对同频干扰进行抑制，最后采用时频二维相关提取目标时延、多普勒频率，从而估计目标DOA(目标到达角度)参数。然而，此种方法在目标信号十分微弱的情况下将失效，同时，由于同频干扰也可能存在于同一个波束以内，这将导致无法提取目标时延、多普勒频率参数，从而无法估计目标的DOA。

发明内容

本发明所要解决的技术问题是，提供一种存在同频干扰情况下有效实现微弱目标回波信号检测的方法。

本发明为解决上述技术问题所采用的技术方案是，一种基于主量分析PCA的外辐射源微弱信号检测方法，包括以下步骤：

步骤1：接收方对接收到的待处理的主通道阵列数字调制信号进行带通采样得到采样序列x(t)，t表示采样时刻；

步骤2：求采样序列x(t)的自相关矩阵R=E{x(t)x^H(t)}，E{·}表示求均值；

步骤3：对自相关矩阵R进行奇异值分解SVD得到R=QΛQ^H，利用分解得到的对角阵Λ以及酉矩阵Q对采样序列x(t)进行线性变换得到接收信号从强到弱的排列y(t)，其中y(t)=[y₁(t) y₂(t) … y_M(t)]^T，M为接收天线阵元个数；

步骤4：将接收信号y(t)中第3路至第M路信号分别与第1路信号y₁(t)做时频分布二维相关得到对应相关值|Ψ_i(τ,f)|，i=3,…,M，L为带通采样长度，τ为时延、f多普勒频移；

步骤5：分别在3至M各路信号对应的相关值|Ψ_i(τ,f)|中搜素峰值点，将各路信号峰值点对应的时延与多普勒频移作为一个有效信号，计算各路信噪比，当有某一路信噪比大于等于门限值时，则将其峰值点对应的时延与多普勒频确定为一个目标的时延与多普勒频移，否则判定该路信号中不存在目标。

本发明的有益效果是，本发明能够实现在强同频干扰以及低信噪比条件下提取微弱目标信号。另外本发明在在时频二维相干匹配提取时延、多普勒频移时，不需要再采用参考天线接收信号作为参考，从而起到了节约硬件成本的作用，另外本发明不用采用自适应对消等步骤，整体运算量也较传统方法有较大的减小。

附图说明

图1：基于广播外辐射源定位的***原理示意图。

图2：传统方法信号处理流程图。

图3：本发明的信号处理流程图。

图4：目标衰减-60dB，干扰衰减-10dB的情况下，传统方法目标信号的时频二维相关图。

图5：目标衰减-60dB，干扰衰减-10dB的情况下，本发明目标信号的时频二维相关图。

图6：目标衰减-60dB，干扰衰减-10dB的情况下，传统方法目标信号的时频二维相关图时延维。

图7：目标衰减-60dB，干扰衰减-10dB的情况下，本发明目标信号的时频二维相关图时延维。

图8：目标衰减-60dB，干扰衰减-10dB的情况下，传统方法目标信号的时频二维相关图频率维。

图9：目标衰减-60dB，干扰衰减-10dB的情况下，本发明目标信号的时频二维相关图频率维。

具体实施方式

本发明处理信号模型如下：

对阵列天线接收到的数字调制信号y(t)进行建模，主要针对带通信号进行处理，设s(t)表示某一个广播电视塔发射的调频广播连续的经过调制后的信号，则为目标的回波信号，其中M为接收天线阵元个数，a₁,…,a_M表示目标回波的衰减，τ₁,…τ_M表示目标的时延，f₁,…,f_M表示目标的多普勒频移，则由阵列天线接收到的信号为 $y\left(t\right)=s\left(t\right)+a_{1}s(t-{\mathrm{\τ}}_1)e^{-j2\mathrm{\π}{f}_{1}t}+\·\·\·+a_{M}s(t-{\mathrm{\τ}}_M)e^{-j2\mathrm{\π}{f}_{M}t}+z\left(t\right)+n\left(t\right),$ z(t)表示同频干扰信号，n(t)表示接收机的热噪声，则本发明所要做的工作是从接收到的信号中提取目标的时延、多普勒频率参数，这是为了在后续处理中能够估计目标DOA，从而最终实现对目标的定位与跟踪做准备。

实施步骤如图3所示：

步骤1：接收方对接收到的待处理的主通道阵列数字调制信号进行带通采样，得到采样序列x(t)=[x₁(t) x₂(t) … x_M(t)]^T，其中M表示阵元数，t表示采样时刻，采样长度为L，T表示矩阵转置。x(t)包含有目标时延、多普勒频率以及DOA的信息。

步骤2：求采样序列x(t)的自相关矩阵R=E{x(t)x^H(t)}，E{·}表示求均值。可假设采样序列在采样长度L内是平稳的随机过程，这样求统计平均可以采用求时间平均代替，则 $R=\frac{1}{L}\underset{t}{\mathrm{\Σ}}x\left(t\right)x^H\left(t\right),$ H表示共轭转置。

步骤3：对自相关矩阵R进行奇异值分解SVD得到R=QΛQ^H，利用分解得到的对角阵Λ以及酉矩阵Q，对x(t)进行线性变换得到其中y(t)=[y₁(t) y₂(t) … y_M(t)]^T，y(t)是将接收信号从强到弱的排列，y₁(t)为最强的直达波信号，即参考信号。

步骤4：采用时频二维相干匹配滤波提取目标的时延和多普勒频移：由于直达波信号比较强，故在进行主分量分析后所得到的第一路信号y₁(t)可作为直达波信号，y₁(t)=Aexp{j[2πf_c(t-D₁)+θ(t-D₁)]}+n(t)，其中n(t)为干扰项，f_c为广播信号的载波频率，D₁为广播信号从发射到接收所经过的时间，θ(t-D₁)为音频信号调制项。而由于干扰功率仅次于直达波功率，故第二路信号y₂(t)可看作是同频干扰信号。y_i(t),i≠1,2为包含目标回波参数的信号，由于经过主分量分析，其中所含干扰和直达波分量已经比较小，故其表达式可写作y_i(t)=Bexp{j[2π(f_c+f_dm)(t-D_m)+θ(t-D_m)]}+m(t)，其中f_dm为目标多普勒频率，D_m为广播信号从发射、经过目标的散射到达雷达接收天线所经过的时间，m(t)可看作与估计所需(τ,f_d)参数无关的干扰项，对y₁(t)做时延可得：

y₁(t-n△t)=Aexp{j[2πf_c(t-D₁-n△t)+θ(t-D₁-n△t)]}+n((t-n△t))

延时长度△t用采样周期，或者采样周期的整数倍。用y₁(t-n△t)的共轭式与y_i(t),i≠1,2相乘可得：

$y_1^{*}(t-\mathrm{n\Δt})y_i\left(t\right)$

$=\mathrm{ABexp}\left\{j\right[2\mathrm{\π}{f}_{\mathrm{dm}}t+2\mathrm{\π}{f}_c(D_1+\mathrm{n\Δt}-D_m)-2\mathrm{\π}{f}_d{D}_m+\mathrm{\θ}(t-D_m)-\mathrm{\θ}(t-D_1-\mathrm{n\Δt})\left]\right\}+N\left(t\right)$

$=\mathrm{ABexp}\left\{j\right[2\mathrm{\π}{f}_{\mathrm{dm}}t+\mathrm{\θ}(t-D_m)-\mathrm{\θ}(t-D_1-\mathrm{n\Δt})+{\mathrm{\θ}}_1\left]\right\}+N\left(t\right)$

*表示共轭矩阵，可以看到，通过调整直达波信号的时延，当D₂=D₁+n△t时，上式中指数项中的两项相位将相消，因而上式将变为多普勒频移项和噪声项，即：

$y_1^{*}(t-\mathrm{n\Δt})y_i\left(t\right)=\mathrm{ABexp}\left\{j\right[2\mathrm{\π}{f}_{\mathrm{dm}}t+{\mathrm{\θ}}_1\left]\right\}+N\left(t\right)$

对上式做FFT运算，即可得到目标的时延和多普勒频移。

将步骤3得到的M路信号的第3路到第M路信号与第1路信号y₁(t)做时频二维相关，i=3,…,M，(τ,f)为第i路信号对应的时延和多普勒频率。当存在目标时，第i路上的|Ψ_i(τ,f)|会产生一个峰值点|Ψ_i(τ_d,f_d)|，搜索到峰值点后，将(τ_d,f_d)作为一个有效信号，计算信噪比， 表示第i路信号上除|Ψ_i(τ_d,f_d)|外，所有|Ψ_i(τ,f)|的均值，若其值大于预先设置的门限值13dB，则表示目标存在，对应的(τ_d,f_d)即为一个目标的时延和多普勒频率；若其值小于这一门限值，则表示第i路信号上不存在目标。

为了验证本发明效果，分别采用两段独立的随机音频信号作为调制信号用以产生调频广播信号，以此作为仿真场景中的直达波信号和同频干扰，载波频率为80kHz，最大频偏为±75kHz，采样频率为200kHz。仿真平台操作***是WIN7，采用MATLAB进行编译。

表1给出混合信号的各仿真参数，并比较在不同条件下传统方法和本发明提出方法的对目标信号的检测性能，最后还比较了传统方法和本发明所提方法的运算量。

表1仿真信号参数

信号	时延（μs）	功率衰减（dB）	到达角（°）	多普勒频率（Hz）
					主台直达波s1	0	0	-50	0
主台多径s2	75	-10	-20	0
					主台多径s3	40	-20	30	0
临台干扰s4	0	-10～-30	8	0
					目标回波s5	439.12	-60～-85	5	119.37

接收阵列阵元数为12。附加的噪声为-50dB，高于目标回波信号10dB～35dB。时频二维相干积累时间为1s。数字波束形成DBF采用切比雪夫加权，波束范围为从-45°到45°，波束间隔为10°。恒虚警检测器采用单元选大恒虚警检测器，设置恒虚警门限系数为4.5。对传统方法性能如表2，本发明提出方法性能如表3：

表2传统方法检测性能表

表3本发明所提方法检测性能表

设目标时延、多普勒频率处时频二维相关最大峰值为A，除通过恒虚警率处理CFAR可检测到的目标峰值以外，时频二维相关图上底噪的平均为B，则上表中最大峰值区分度定义为10log(A/B)。

通过表2和表4所提供的仿真数据可以说明本发明所提出的方法在有强同频干扰，目标信号和干扰在同一个波束以内，以及目标信号十分微弱的情况下都能有效的检测到目标信号，并提取出目标的时延、多普勒频率参数。

目标衰减-60dB，干扰衰减-10dB的条件下，传统方法和本发明的目标信号时频二维相关图分别为图4、图5；传统方法和本发明的目标信号时频二维相关图时延维分别为图6、图7；传统方法和本发明的目标信号时频二维相关图频率维分别为图8、图9。可以分别从，传统方法和本发明的时频二维相关图时延维和多普勒维观察对比目标，通过对比可以发现，传统方法时频二维相关图的底噪很高，目标的最大峰值区分度较低，而本发明时频二维相关图的底噪较小，且目标最大峰值区分度较高，检测效果好。

由于传统方法和本发明都要用到时频二维相关，故比较计算量只需比较时频二维相关前的运算量，设接收阵列阵元个数为M，快拍数为N，波束个数为L，则比较结果如表4、表5所示：

表4传统方法运算量表

表5本发明所提方法运算量表

通过以上计算对比可以看出，本发明比传统方法的运算量有一定量的减小。且本方法不需要使用参考天线对准直达波信号，在一定程度上节省了硬件的成本。

Claims (5)

1.一种基于PCA的外辐射源微弱信号检测方法，包括以下步骤：

步骤1：接收方对接收到的待处理的主通道阵列数字调制信号进行带通采样得到采样序列x(t)，t表示采样时刻；

步骤2：求采样序列x(t)的自相关矩阵R；

步骤3：对自相关矩阵R进行奇异值分解SVD得到R＝QΛQ^H，利用分解得到的对角阵Λ以及酉矩阵Q对采样序列x(t)进行线性变换得到接收信号从强到弱的排列y(t)，其中y(t)＝[y₁(t)  y₂(t)…y_M(t)]^T，M为接收天线阵元个数；

步骤4：将接收信号y(t)中第3路至第M路信号分别与第1路信号y₁(t)做时频分布二维相关得到对应相关值|Ψ_i(τ,f)|， $\left|{\mathrm{\Ψ}}_i(\mathrm{\τ},f)\right|=\left|\underset{t=0}{\overset{L-1}{\mathrm{\Σ}}}{y}_i\left(t\right){y_1}^{*}(t-\mathrm{\τ})e^{j2\mathrm{\πf}/L}\right|,i=3,...,M,$ L为带通采样长度，τ为时延、f多普勒频移；

步骤5：分别在3至M各路信号对应的二维相关值|Ψ_i(τ,f)|中搜素峰值点，将各路信号峰值点对应的时延与多普勒频移作为一个有效信号，计算各路信噪比，当有某一路信噪比大于等于门限值时，则将其峰值点对应的时延与多普勒频移确定为一个目标的时延与多普勒频移，否则判定该路信号中不存在目标。

2.如权利要求1所述一种基于PCA的外辐射源微弱信号检测方法，其特征在于，自相关矩阵R＝E{x(t)x^H(t)}，E{·}表示求均值。

3.如权利要求1所述一种基于PCA的外辐射源微弱信号检测方法，其特征在于， $R=\frac{1}{L}\underset{t}{\mathrm{\Σ}}x\left(t\right)x^H\left(t\right),$ H表示共轭转置。

4.如权利要求1、2或3所述一种基于PCA的外辐射源微弱信号检测方法，其特征在于，步骤5中计算各路信噪比的方法为：|Ψ_i(τ_d,f_d)|为第i路上的二维相关值|Ψ_i(τ,f)|的峰值点，表示第i路信号上除|Ψ_i(τ_d,f_d)|外，所有|Ψ_i(τ,f)|的均值。

5.如权利要求4所述一种基于PCA的外辐射源微弱信号检测方法，其特征在于，所述门限值为13dB。