GNSS-R外辐射源雷达微弱目标回波盲检测方法

温媛媛; 白琳; 尚社; 宋大伟; 王俊; WEN Yuanyuan; BAI Lin; SHANG She; SONG Dawei; WANG Jun

doi: 10.11887/j.cn.202405013

温媛媛¹ ，白琳² ，尚社¹ ，宋大伟¹ ，王俊³

1. 中国空间技术研究院西安分院空间微波通信全国重点实验室, 陕西西安 710100

2. 中国空间技术研究院西安分院, 陕西西安 710100

3. 西安电子科技大学雷达信号处理国家重点实验室, 陕西西安 710071

基金项目: 重点实验室基金资助项目（6142411332103，HTKJ2022KL504005）

详细信息

作者简介

温媛媛（1980—），女，河南洛阳人，高级工程师，博士，E-mail: wyy031980@163.com

中图分类号: TN958

文献标识码: A

文章编号: 1001-2486(2024)05-121-10

Blind detection method of weak target echo in GNSS-R passive radar

WEN Yuanyuan¹ ， BAI Lin² ， SHANG She¹ ， SONG Dawei¹ ， WANG Jun³

1. National Key Laboratory of Science and Technology on Space Microwave, China Academy of Space Technology (Xi′an), Xi′an 710100 , China

2. China Academy of Space Technology (Xi′an), Xi′an 710100 , China

3. National Key Laboratory of Radar Signal Processing, Xidian University, Xi′an 710071 , China

摘要

针对传统全球导航卫星系统反射(global navigation satellite system-reflected，GNSS-R)外辐射源雷达系统所采用的双通道配置结构在信号处理上运算量较大且在工程实现上硬件成本较高的问题, 提出一种基于高阶循环累积量的单通道GNSS-R外辐射源雷达微弱目标回波盲检测方法。利用主成分分析方法将强直达波信号从单通道接收的混合信号中抽取出来; 通过直达波信号的高阶循环频率估计出微弱回波信号的高阶循环频率; 基于微弱回波信号高阶循环频率与其他信号不同的特性将其从单通道混合信号中提取出来，通过匹配滤波实现目标检测。实验结果表明, 所提方法在没有先验信息的情况下能有效提取出微弱目标的回波信号, 且与传统二维匹配滤波目标检测方法相比有更好的目标检测性能。

关键词

外辐射源雷达 / 双通道配置 / 高阶循环累积量 / 单通道 / 主成分分析

Abstract

Aiming at the problem that the dual-channel configuration of the traditional GNSS-R(global navigation satellite system-reflected) passive radar system has a large amount of computation in signal processing and a high hardware cost in engineering implementation, a blind detection method of the weak target echo of single-channel GNSS-R passive radar based on higher-order cyclic cumulant was presented. Strong direct wave signal was extracted from the single-channel mixed signal using the principal component analysis method. Higher-order cyclic frequency of weak target echo was estimated by higher-order cyclic frequency of direct wave signal. Weak echo signal was extracted from the single-channel mixed signal based on its different characteristics of higher-order cyclic frequency from other signals, so that the object detection was realized by matched filtering. Simulation results show that the proposed method can effectively extract the weak target echo without prior information, and has better target detection performance compared with the traditional two-dimensional matched filter target detection method.

Keywords

passive radar / dual-channel configuration / higher-order cyclic cumulant / single-channel / principal component analysis

1 单通道GNSS-R外辐射源雷达微弱目标回波盲检测模型 2 单通道GNSS-R外辐射源雷达微弱目标回波盲检测方法 2.1 基于PCA的直达波抽取 2.2 高阶循环频率估计 2.3 基于高阶循环平稳度的微弱目标盲提取 3 仿真分析 3.1 实验一：算法盲提取性能分析 3.2 实验二：算法目标检测性能分析 4 结论

全球导航卫星系统（global navigation satellite system，GNSS）利用美国的GPS、俄罗斯的格洛纳斯（GLONASS）、欧洲正在研发中的伽利略（GALILEO）和中国的北斗（BD）等卫星导航系统中的一个或多个系统进行导航定位，并同时提供卫星的完备性检验信息和足够的导航安全性告警信息^[1-3]。GNSS具备全天候连续提供全球高精度导航的能力，具有发射频率稳定、授时精确易于同步、直达波恢复容易等特点^[4]。针对GNSS卫星信号的优势，利用GNSS作为外辐射源可以有效解决双基雷达探测系统中的同步问题。全球导航卫星系统反射（global navigation satellite system-reflected，GNSS-R）雷达属于外辐射源雷达体制，可利用GNSS-R信号作为非合作的照射源信号进行无源探测。

在传统的主动雷达探测系统中，系统本身具有发射机，则参考信号已知。但在GNSS-R外辐射源雷达系统中，参考信号无法直接获得，且阵列天线接收到来自目标反射的回波信号能量远低于来自导航卫星的强直达波和经地面及建筑物反射的多径杂波与噪声信号，很难实现目标的探测^[5]。因此要想准确提取出微弱的目标回波进一步实现目标探测，一般情况下必须进行双通道配置。参考通道用指向信号源发射方向的天线接收外辐射源的直达波作为参考信号；回波通道用阵列天线接收混合信号作为回波信号，其中混合信号包括经过目标所反射回来的信号、外辐射源的直达波信号及多径和噪声等起到扰乱作用的信号^[6-7]。通过对参考通道和回波通道里的信号分别进行处理，可以实现对目标的探测。一般情况下，参考通道需要进行直达波提纯，目前使用比较多的方法是盲均衡方法和信号重构方法。文献^[8-11]研究了外辐射源雷达系统中参考通道直达波的提纯方法。回波通道一般需要进行干扰信号抑制，目前较为常用的方法是自适应滤波算法。文献^[12-14]研究了外辐射源雷达系统中回波通道一些干扰抑制方法。传统GNSS-R外辐射源雷达系统双通道配置：一方面信号处理过程较为复杂，运算量相对较大；另一方面双通道接收通道结构较为复杂，硬件成本较高。在实际场合中，用指向GNSS的天线接收直达波信号作为参考信号时，需保证GNSS方向完全已知。但在电磁环境日益复杂的今天，尤其是复杂的战场环境中，GNSS-R外辐射源雷达接收机周围可能存在除GNSS之外的其他外辐射源，这必然导致外辐射源基本信息未知，因此直接用上述的指向性极强的天线接收GNSS信号作为参考信号的方案存在缺陷。

针对上述问题，本文利用GNSS-R外辐射源雷达系统中信号的高阶循环平稳特性，提出一种单通道GNSS-R雷达系统微弱目标回波盲检测方法。该方法包含三个步骤：首先利用主成分分析（principal component analysis，PCA）方法将强直达波信号从回波通道接收信号中提取出来，将提取出的直达波信号作为参考信号；然后利用信号的高阶循环累积量特性估计出GNSS直达波信号、多径及目标回波信号高阶循环频率；最后将直达波信号、微弱回波信号及多径信号的高阶循环频率的不同作为先验信息，利用盲提取方法将微弱回波信号从单通道混合信号中提取出来，进一步将微弱回波信号与参考信号进行匹配滤波，得到目标的检测结果。

1 单通道GNSS-R外辐射源雷达微弱目标回波盲检测模型

在传统GNSS-R外辐射源雷达系统中，参考通道用指向性极强的天线接收直达波信号。从经济角度来说，多增加一个通道就需要增加天线和接收机，额外增加了硬件成本。参考通道中直达波信号的纯净度对后续的微弱目标探测至关重要，因此传统GNSS-R外辐射源雷达系统进行微弱目标探测时首先得对参考通道接收到的信号进行干扰抑制、杂波相消及直达波提纯等相关的处理，这些处理过程也都会额外增加运算量。

本文提出的单通道GNSS-R外辐射源雷达系统打破传统外辐射源雷达的双通道信号接收模型，采用单通道形式接收信号，去掉参考通道就可以减少一个天线和一个接收机，节约了硬件成本，后续信号处理过程只对单通道接收到的信号进行处理即可。单通道所接收到的信号分量包括直达波信号、多径杂波信号、目标回波信号、噪声信号及其他未知干扰信号，单通道GNSS-R外辐射源雷达微弱目标回波盲检测模型如图1和图2所示。

图1GNSS-R外辐射源雷达微弱目标盲检测模型

Fig.1GNSS-R passive radar weak target echo blind detection model

在单通道GNSS-R外辐射源雷达系统中，假设接收天线有M路接收通道，单通道中接收到的信号数学模型如式（1）所示。

图2单通道微弱目标回波盲检测信号处理流程

Fig.2Signal processing flow chart of single channel weak target echo blind detection

X (t) = H \cdot \sum_{i = 1}^{N} r (θ_{i}) s_{i} (t) + n (t) = A \cdot S (t) + n (t)

(1)

式中：X（t）为回波通道天线接收到的观测信号，即在处理器前的数据；A为M×N混叠矩阵，它为天线阵的阵列流形，即天线阵的响应函数[r（θ₁），···，r（θ_N）]和接收天线阵列阵元耦合矩阵H的乘积；

S （ t ） = [A_{0} s （ t ） ， \dots ， b_{i} s (t - τ_{i}) e^{j 2 {π f}_{i} t} ， \dots ， {A_{1} s (t - τ_{t r}) e^{j 2 π f_{d r} t}]}^{T}

为源信号（包括直达波信号、多径干扰及回波信号），其中A₀为直达波信号的幅度，A₁=ξA₀为目标回波信号的幅度，ξ为目标回波的能量散射系数，τ_tr为目标回波的时延，f_dr为目标回波的多普勒频率，b_i为第i（1≤i≤P）个多径杂波信号的幅度，τ_i为第i个多径杂波的时延，f_i为第i个多径杂波的多普勒频率，P为多径干扰个数；n（t）为M×1维噪声信号，包括外部噪声、内部热噪声及电噪声等，通常视为高斯白噪声。在下面算法中，假设各信号间相互独立且与噪声亦相互独立。

2 单通道GNSS-R外辐射源雷达微弱目标回波盲检测方法

2.1 基于PCA的直达波抽取

PCA是一种应用最为广泛的统计数据分析、特征提取和数据压缩的经典方法。其基本思想是提取空间原始数据中的主要特征，减少数据冗余，寻找原始变量冗余度最小的一个子集来更好地表示原始变量。具体方法是通过在样本空间中找到一组正交向量，即主成分，进而最大限度地表示原始数据中的信息，从而实现对多维矢量数据的主要特征分量提取^[15-16]。PCA通过对原始数据自相关矩阵进行处理，将自相关矩阵奇异值分解后变换为一组按照奇异值大小进行降序排列的主成分，各主成分之间互不相关，再提取出低阶主成分进行信号重构以去除不相关的噪声，提高信号的信噪比^[17]。

单通道GNSS-R外辐射源雷达系统接收通道中信号表示为x_i（t_j），i=1，2，···，M，t_j为采样时间，j=1，2，···，p，其自相关矩阵Γ的矩阵元素γ_jk表示为：

γ_{j k} = \sum_{i = 1}^{M} x_{i} (t_{j}) x_{i} (t_{k}) j, k = 1,2, \dots, p

(2)

利用奇异值分解法，将自相关矩阵Γ分解为特征值与特征向量的形式：

Γ = R Λ R^{T}

(3)

式中：Λ为由大到小排列的特征值对角矩阵； R=[r_jk]_j_，_k_{=1，2，···，}_p为各个特征值对应的特征向量矩阵，代表数据能量递减的各个主成分分量的方向，R^T又称为旋转矩阵。

单通道GNSS-R外辐射源雷达系统接收通道数据利用旋转矩阵R^T进行线性映射，得到能量递减的各个主成分方向的系数，称为第1，2，···，p个主成分：

Ψ = R^{T} X 或 ψ_{k i} = \sum_{j = 1}^{p} r_{j k} x_{i} (t_{j})

(4)

式中，Ψ为主成分矩阵，ψ_ki为第i组GNSS-R雷达接收数据的第k个主成分。

根据矩阵计算方法，GNSS-R外辐射源雷达接收数据X可由主成分Ψ与特征向量矩阵R进行重构，即：

X = R Ψ 或 x_{i} (t_{j}) = \sum_{j = 1}^{p} r_{j k} ψ_{k i}

(5)

由于GNSS-R外辐射源雷达接收数据的能量主要集中在低阶主成分上，而高阶主成分代表着不相关噪声，因此通过少量低阶（L个）主成分重构GNSS-R雷达接收数据（见式（6）），既能保证主成分信号（直达波）精度，又能达到去除不相关噪声的目的。

{\tilde{x}}_{i} (t_{j}) = \sum_{k = 1}^{L} r_{j k} ψ_{k i}

(6)

式中，

{\tilde{x}}_{i} (t_{j})

为采用少量低阶主成分进行重构的GNSS-R外辐射源雷达系统直达波信号。

2.2 高阶循环频率估计

单回波通道中接收到的混合信号包含时频域都相互重叠的直达波、微弱目标回波及多径干扰，很难将微弱的目标回波从混合信号中提取出来。目标回波信号和其他信号的差异只体现在高阶循环频率上，因此，考虑利用高阶循环频率来抽取目标回波信号。通过计算信号的高阶循环累积量可得到信号的高阶循环频率，高阶循环累积量理论上可以完全抑制任何平稳（高斯或非高斯）噪声，故本文方法利用三阶循环累积量函数进行循环频率的估计。在单通道GNSS-R外辐射源雷达系统中，一个接收天线中的信号可表示为：

\begin{matrix} x (t) = A_{0} s (t) + \sum_{i = 1}^{P} b_{i} s (t - τ_{i}) e^{j 2 π f_{i} t} + \\ A_{1} s (t - τ_{tr}) e^{j 2 π f_{d r} t} + n (t) \end{matrix}

(7)

系统单通道接收信号x（t）的三阶循环累积量函数可表示为：

C_{3 x}^{β} (τ_{1}, τ_{2}) = {⟨x (t) x (t + τ_{1}) x (t + τ_{2}) e^{- j 2 π β t}⟩}_{t}

(8)

式中，β为三阶循环频率，τ₁、τ₂为两个时延，〈·〉_t表示时间平均。三阶循环累积量是一个三维矩阵，其提供信号可分离性的这个性质不会因为时间延迟量的不同而改变，时间延迟量的不同也不会影响分离信号算法的性能，于是本文只取时间延迟都为0的方向上的一个值作为代表。将式（7）代入式（8）进行计算，当τ₁=τ₂=0时有：

\begin{matrix} C_{3 x}^{β} (0) = {⟨x (t) x (t) x (t) e^{- j 2 π β t}⟩}_{t} = A_{0}^{3} {⟨s^{3} (t) e^{- j 2 π β t}⟩}_{t} + \\ 3 \sum_{i = 1}^{P} A_{0}^{2} b_{i} {⟨s^{2} (t) s (t - τ_{i}) e^{j 2 π f_{i} t} e^{- j 2 π β t}⟩}_{t} + \\ 3 A_{0}^{2} A_{1} {⟨s^{2} (t) s (t - τ_{tr}) e^{- j 2 π (β - f_{d r}) t}⟩}_{t} + \\ 3 \sum_{i = 1}^{P} \sum_{j = 1}^{P} A_{0} b_{i} b_{j} {⟨s (t) s (t - τ_{i}) s (t - τ_{j}) e^{j 2 π f_{i} t} e^{j 2 π f_{j} t} e^{- j 2 π β t}⟩}_{t} + \\ 6 \sum_{i = 1}^{P} A_{0} A_{1} b_{i} {⟨s (t) s (t - τ_{i}) s (t - τ_{tr}) e^{j 2 π f_{i} t} e^{- j 2 π (β - f_{d r}) t}⟩}_{t} + \\ 3 A_{0} A_{1}^{2} {⟨s (t) s^{2} (t - τ_{tr}) e^{- j 2 π (β - 2 f_{d r}) t}⟩}_{t} + \end{matrix}

\begin{matrix} \sum_{i = 1}^{P} \sum_{j = 1}^{P} \sum_{k = 1}^{P} b_{i} b_{j} b_{k} {⟨s (t - τ_{i}) s (t - τ_{j}) s (t - τ_{k}) e^{j 2 π f_{i} t} e^{j 2 π f_{j} t} e^{j 2 π f_{k} t} e^{- j 2 π β t}⟩}_{t} + \\ 3 \sum_{i = 1}^{P} \sum_{j = 1}^{P} A_{1} b_{i} b_{j} {⟨s (t - τ_{i}) s (t - τ_{j}) s (t - τ_{t r}) e^{j 2 π f_{i} t} e^{j 2 π f_{j} t} e^{- j 2 π (β - f_{d r}) t}⟩}_{t} + \\ 3 \sum_{i = 1}^{P} A_{1}^{2} b_{i} {⟨s (t - τ_{i}) s^{2} (t - τ_{tr}) e^{j 2 π f_{i} t} e^{- j 2 π (β - 2 f_{d r}) t}⟩}_{t} + \\ A_{1}^{3} {⟨s^{3} (t - τ_{tr}) e^{- j 2 π (β - 3 f_{d r}) t}⟩}_{t} \\ = A_{0}^{3} C_{3 s}^{β} (0, 0) + 3 \sum_{i = 1}^{P} A_{0}^{2} b_{i} C_{3 s}^{β - f_{i}} (0, - τ_{i}) + \end{matrix}

\begin{matrix} 3 A_{0}^{2} A_{1} C_{3 s}^{β - f_{d r}} (0, - τ_{d r}) + 3 \sum_{i = 1}^{P} \sum_{j = 1}^{P} A_{0} b_{i} b_{j} C_{3 s}^{β - f_{i} - f_{j}} (- τ_{i}, - τ_{j}) + \\ 6 \sum_{i = 1}^{P} A_{0} A_{1} b_{i} C_{3 s}^{β - f_{d r} - f_{i}} (- τ_{i}, - τ_{t r}) + 3 A_{0} A_{1}^{2} C_{3 s}^{β - 2 f_{d r}} (- τ_{t r}, - τ_{t r}) + \\ \sum_{i = 1}^{P} \sum_{j = 1}^{P} \sum_{k = 1}^{P} b_{i} b_{j} b_{k} C_{3 s}^{β - f_{i} - f_{j} - f_{k}} (τ_{i} - τ_{j}, τ_{i} - τ_{k}) + \\ 3 \sum_{i = 1}^{P} \sum_{j = 1}^{P} A_{1} b_{i} b_{j} C_{3 s}^{β - f_{d r} - f_{i} - f_{j}} (τ_{i} - τ_{j}, τ_{i} - τ_{t r}) + \\ 3 \sum_{i = 1}^{P} A_{1}^{2} b_{i} C_{3 s}^{β - 2 f_{d r} - f_{i}} (0, τ_{t r} - τ_{i}) + A_{1}^{3} C_{3 s}^{β - 3 f_{d r}} (0,0) \end{matrix}

(9)

由式（9）可见，函数峰值出现在频率值为0、α、α+f_dr、α+2f_dr、α+3f_dr、f_dr、2f_dr、3f_dr处，其中α为直达波的三阶循环频率。由此可得到直达波的三阶循环频率α与目标回波的多普勒频移f_dr。式（9）结果中的最后一项

A_{1}^{3} C_{3 s}^{β - 3 f_{d r}} （ 0，0 ）

表示目标回波信号的三阶循环累积量函数，因此目标回波信号的三阶循环频率为α+3f_dr、3f_dr。

2.3 基于高阶循环平稳度的微弱目标盲提取

1991年，Živanović等^[18]提出了循环平稳度（degree of cyclostationarity，DCS）的概念，用信号x在循环频率α处的循环自相关函数的积分与其自相关函数积分的比值来度量信号在其循环频率α处的循环平稳程度。DCS的值界于[0，1]，当α≠0的时候，DCS的值越大则表明信号在该循环频率处的周期平稳程度越强。本文将这个循环平稳度的概念扩展到高阶循环平稳信号的应用中，提出三阶循环平稳度的定义。DCS是用来度量信号在循环频率α处的循环平稳程度。连续信号的三阶循环平稳度定义^[19-20]为：

D C S^{3 α} = \frac{\int_{- \infty}^{+ \infty} \int_{- \infty}^{+ \infty} {|C_{3 x}^{α} (τ_{1}, τ_{2})|}^{2} d τ_{1} d τ_{2}}{\int_{- \infty}^{+ \infty} \int_{- \infty}^{+ \infty} {|C_{3 x}^{0} (τ_{1}, τ_{2})|}^{2} d τ_{1} d τ_{2}}

(10)

取该函数作为本文单通道GNSS-R雷达系统微弱回波信号盲提取的准则。在盲提取过程中计算DCS³^α的最大值，就能够较好地提取出混合信号中的微弱回波信号。

为了验证本文盲提取准则的有效性，考虑单通道GNSS-R外辐射源雷达系统接收信号中混合有三个源（s₁（t）、s₂（t）和s₃（t））的情况，混合模型假定是线性的，微弱回波信号的估计量y₁（t）是三个源的组合，如下式所示：

y_{1} (t) = a (θ) s_{1} (t) + b (θ) s_{2} (t) + c (θ) s_{3} (t)

(11)

式中系数a（θ）、b（θ）及c（θ）的大小最终确定出旋转矩阵，旋转矩阵对混合信号进行旋转择优，最后得到微弱回波信号。

计算微弱回波信号的估计量y₁（t）在其三阶循环频率α处的三阶循环累积量：

\begin{matrix} C_{3 y_{1}}^{α} (τ_{1}, τ_{2}) = {⟨y_{1} (t) y_{1} (t + τ_{1}) y_{1} (t + τ_{2}) e^{- j 2 π α t}⟩}_{t} \\ = {⟨(Ω_{1} + Ω_{2} + Ω_{3}) e^{- j 2 π α t}⟩}_{t} \\ = a^{3} (θ) C_{3 s_{1}}^{α} (τ_{1}, τ_{2}) + b^{3} (θ) C_{3 s_{2}}^{α} (τ_{1}, τ_{2}) + \\ c^{3} (θ) C_{3 s_{3}}^{α} (τ_{1}, τ_{2}) \end{matrix}

(12)

式中，

Ω_{1} = a^{3} (θ) s_{1} (t) s_{1} (t + τ_{1}) s_{1} (t + τ_{2})

\begin{matrix} Ω_{2} = b^{3} (θ) s_{2} (t) s_{2} (t + τ_{1}) s_{2} (t + τ_{2}) \\ Ω_{3} = c^{3} (θ) s_{3} (t) s_{3} (t + τ_{1}) s_{3} (t + τ_{2}) \end{matrix}

如果本文算法能很好地实现微弱回波信号的盲抽取，即y₁（t）是源信号s₁（t）的估计量，则y₁（t）的三阶循环频率α也是s₁（t）的三阶循环频率，在此三阶循环频率处信号s₂（t）和s₃（t）的三阶循环累积量为0。故由式（12）最后可推导出以下结果：

C_{3 y_{1}}^{α} (τ_{1}, τ_{2}) = a^{3} (θ) C_{3 s_{1}}^{α} (τ_{1}, τ_{2})

(13)

将α=0代入式（12）可得：

\begin{matrix} C_{3 y_{1}}^{0} (τ_{1}, τ_{2}) = a^{3} (θ) C_{3 s_{1}}^{0} (τ_{1}, τ_{2}) + b^{3} (θ) C_{3 s_{2}}^{0} (τ_{1}, τ_{2}) + \\ c^{3} (θ) C_{3 s_{3}}^{0} (τ_{1}, τ_{2}) \end{matrix}

(14)

将式（13）和式（14）代入式（10），可得：

\begin{matrix} f (θ) = D C S^{3 α} = \frac{\int_{- \infty}^{+ \infty} \int_{- \infty}^{+ \infty} {|a^{3} (θ) C_{3 s_{1}}^{α} (τ_{1}, τ_{2})|}^{2} d τ_{1} d τ_{2}}{\int_{- \infty}^{+ \infty} \int_{- \infty}^{+ \infty} Π d τ_{1} d τ_{2}} \\ = \frac{\int_{- \infty}^{+ \infty} \int_{- \infty}^{+ \infty} a^{6} (θ) {|C_{3 s_{1}}^{α} (τ_{1}, τ_{2})|}^{2} d τ_{1} d τ_{2}}{\int_{- \infty}^{+ \infty} \int_{- \infty}^{+ \infty} Π d τ_{1} d τ_{2}} \end{matrix}

(15)

\begin{matrix} \int_{- \infty}^{+ \infty} \int_{- \infty}^{+ \infty} Π d τ_{1} d τ_{2} \\ = \int_{- \infty}^{+ \infty} \int_{- \infty}^{+ \infty} Γ d τ_{1} d τ_{2} + \\ \int_{- \infty}^{+ \infty} \int_{- \infty}^{+ \infty} 2 a^{3} (θ) b^{3} (θ) |C_{3 s_{1}}^{0} (τ_{1}, τ_{2}) C_{3 s_{2}}^{0} (τ_{1}, τ_{2})| d τ_{1} d τ_{2} + \\ \int_{- \infty}^{+ \infty} \int_{- \infty}^{+ \infty} 2 a^{3} (θ) c^{3} (θ) |C_{3 s_{1}}^{0} (τ_{1}, τ_{2}) C_{3 s_{3}}^{0} (τ_{1}, τ_{2})| d τ_{1} d τ_{2} + \\ \int_{- \infty}^{+ \infty} \int_{- \infty}^{+ \infty} 2 b^{3} (θ) c^{3} (θ) |C_{3 s_{2}}^{0} (τ_{1}, τ_{2}) C_{3 s_{3}}^{0} (τ_{1}, τ_{2})| d τ_{1} d τ_{2} \end{matrix}

(16)

其中，

\begin{matrix} Π = ∣ a^{3} (θ) C_{3 s_{1}}^{0} (τ_{1}, τ_{2}) + b^{3} (θ) C_{3 s_{2}}^{0} (τ_{1}, τ_{2}) + \\ {c^{3} (θ) C_{3 s_{3}}^{0} (τ_{1}, τ_{2})|}^{2} \\ Γ = \\ a^{6} (θ) {|C_{3 s_{1}}^{0} (τ_{1}, τ_{2})|}^{2} + b^{6} (θ) {|C_{3 s_{2}}^{0} (τ_{1}, τ_{2})|}^{2} + \\ c^{6} (θ) {|C_{3 s_{3}}^{0} (τ_{1}, τ_{2})|}^{2} \end{matrix}

因为f（θ）是关于θ的函数，与θ无关的部分可以看作常数，所以可以分别令

\int_{- \infty}^{+ \infty} \int_{- \infty}^{+ \infty} {|C_{3 s_{1}}^{α} (τ_{1}, τ_{2})|}^{2} d τ_{1} d τ_{2} = ξ_{1}

(17)

\int_{- \infty}^{+ \infty} \int_{- \infty}^{+ \infty} {|C_{3 s_{1}}^{0} (τ_{1}, τ_{2})|}^{2} d τ_{1} d τ_{2} = ξ_{2}

(18)

\int_{- \infty}^{+ \infty} \int_{- \infty}^{+ \infty} {|C_{3 s_{2}}^{0} (τ_{1}, τ_{2})|}^{2} d τ_{1} d τ_{2} = ξ_{3}

(19)

\int_{- \infty}^{+ \infty} \int_{- \infty}^{+ \infty} {|C_{3 s_{3}}^{0} (τ_{1}, τ_{2})|}^{2} d τ_{1} d τ_{2} = ξ_{4}

(20)

由式（18）和式（19）可得出：

\int_{- \infty}^{+ \infty} \int_{- \infty}^{+ \infty} |C_{3 s_{1}}^{0} (τ_{1}, τ_{2}) C_{3 s_{2}}^{0} (τ_{1}, τ_{2})| d τ_{1} d τ_{2} = \sqrt{ξ_{2} ξ_{3}}

(21)

由式（18）和式（20）可得出：

\int_{- \infty}^{+ \infty} \int_{- \infty}^{+ \infty} |C_{3 s_{1}}^{0} (τ_{1}, τ_{2}) C_{3 s_{3}}^{0} (τ_{1}, τ_{2})| d τ_{1} d τ_{2} = \sqrt{ξ_{2} ξ_{4}}

(22)

由式（19）和式（20）可得出：

\int_{- \infty}^{+ \infty} \int_{- \infty}^{+ \infty} |C_{3 s_{2}}^{0} (τ_{1}, τ_{2}) C_{3 s_{3}}^{0} (τ_{1}, τ_{2})| d τ_{1} d τ_{2} = \sqrt{ξ_{3} ξ_{4}}

(23)

将上面公式代入式（15），可得到f（θ）的简化形式：

\begin{matrix} f (θ) = \frac{a^{6} (θ) ξ_{1}}{Λ + Z} \\ = \frac{a^{6} (θ) ξ_{1}}{{[a^{3} (θ) \sqrt{ξ_{2}} + b^{3} (θ) \sqrt{ξ_{3}} + c^{3} (θ) \sqrt{ξ_{4}}]}^{2}} \end{matrix}

(24)

其中，

\begin{matrix} Λ = a^{6} (θ) ξ_{2} + b^{6} (θ) ξ_{3} + c^{6} (θ) ξ_{4} + 2 a^{3} (θ) b^{3} (θ) \sqrt{ξ_{2} ξ_{3}} \\ Z = 2 a^{3} (θ) c^{3} (θ) \sqrt{ξ_{2} ξ_{4}} + 2 b^{3} (θ) c^{3} (θ) \sqrt{ξ_{3} ξ_{4}} \end{matrix}

只需要f（θ）最大化，即目标信号在α处的三阶循环平稳度最大，即可实现目标信号的盲提取。下面来证明这个准则：

\begin{matrix} \frac{d f (θ)}{d θ} = d \{\frac{a^{6} (θ) ξ_{1}}{{[a^{3} (θ) \sqrt{ξ_{2}} + b^{3} (θ) \sqrt{ξ_{3}} + c^{3} (θ) \sqrt{ξ_{4}}]}^{2}}\} / d θ \\ = \frac{6 a^{5} (θ) ξ_{1} {[a^{3} (θ) \sqrt{ξ_{2}} + b^{3} (θ) \sqrt{ξ_{3}} + c^{3} (θ) \sqrt{ξ_{4}}]}^{2} \frac{d a (θ)}{d θ}}{{[a^{3} (θ) \sqrt{ξ_{2}} + b^{3} (θ) \sqrt{ξ_{3}} + c^{3} (θ) \sqrt{ξ_{4}}]}^{4}} - \\ \frac{6 a^{6} (θ) ξ_{1} [a^{3} (θ) \sqrt{ξ_{2}} + b^{3} (θ) \sqrt{ξ_{3}} + c^{3} (θ) \sqrt{ξ_{4}}] \cdot Φ}{{[a^{3} (θ) \sqrt{ξ_{2}} + b^{3} (θ) \sqrt{ξ_{3}} + c^{3} (θ) \sqrt{ξ_{4}}]}^{4}} \end{matrix}

(25)

其中，

Φ = a^{2} (θ) \sqrt{ξ_{2}} \frac{d a (θ)}{d θ} + b^{2} (θ) \sqrt{ξ_{3}} \frac{d b (θ)}{d θ} + c^{2} (θ) \sqrt{ξ_{4}} \frac{d c (θ)}{d θ}

由上式可看出，当θ=θ_opt时，f（θ）达到最佳分离条件，则应满足b（θ）=0及c（θ）=0，此时有：

{\frac{d f (θ)}{d θ}|}_{θ = θ_{o p t}} = 0

(26)

由此得证将DCS³^α最大化作为微弱回波信号的盲提取准则比较合适。

3 仿真分析

本节利用MATLAB软件对本文方法的盲提取性能及目标检测性能进行相应的仿真分析，并与其他方法进行对比。使用如下所示四个源信号：

\{\begin{matrix} s_{1} (t) = s i n (2 π f_{c} t) + n (t) \\ s_{2} (t) = \{s i n [2 π f_{c} (t - τ_{1})] + n (t - τ_{1})\} e^{j 2 π f_{1} t} \\ s_{3} (t) = \{s i n [2 π f_{c} (t - τ_{2})] + n (t - τ_{2})\} e^{j 2 π f_{2} t} \\ s_{4} (t) = \{s i n [2 π f_{c} (t - τ_{t r})] + n (t - τ_{t r})\} e^{j 2 π f_{d r} t} \end{matrix}

(27)

其中：s₁（t）、s₂（t）、s₃（t）、s₄（t）分别为直达波信号、多径信号1、多径信号2及回波信号（即目标信号y（t））；载波频率f_c=10 MHz； τ₁、τ₂、τ_tr分别为多径信号1、多径信号2及回波信号的时延；f₁、 f₂、 f_dr分别为多径信号1、多径信号2及回波信号的多普勒频率；n（t）为白噪声信号。

3.1 实验一：算法盲提取性能分析

使用相似系数这一常用的盲提取评价指标来定量分析本文算法的提取效果，该系数本质上是相关系数，它的值越大表示算法的盲提取效果越好。相似系数是提取出的信号

\hat{y} （ t ）

与目标信号y（t）相似性的度量，定义为：

ζ (y, \hat{y}) = \frac{|\sum_{t = 1}^{N} \hat{y} (t) y (t)|}{\sqrt{\sum_{t = 1} {\hat{y}}^{2} (t) \sum_{t = 1}^{N} y^{2} (t)}}

(28)

当信干比（signal to interference ratio，SIR）为0 dB时，源信号和用本文方法抽取出的直达波信号及提取出的回波信号波形对比如图3所示，其中图3（a）为源信号，由上至下分别为s₁（t）、s₂（t）、s₃（t）、s₄（t），图3（b）为本文方法恢复出的直达波信号，图3（c）为本文方法恢复出的回波信号，其中纵轴表示幅值，横轴表示采样点。

从图3可看出，本文方法恢复出的回波信号的波形与源信号s₄（t）的波形十分相似，恢复出的直达波信号的波形与源信号s₁（t）的波形也十分相似。为了量化本文方法的性能，在不同SIR条件下，计算提取出的回波信号与源信号s₄（t）的相似系数，结果如表1所示。

图3源信号及用本文方法恢复出的直达波信号和回波信号

Fig.3Source signals and the direct wave signal and echo signal recovered by the proposed method

表1相似系数

Tab.1 Similarity coefficient

从表1可以看出，在SIR较低的情况下，利用本文方法也可以提取出目标信号。接下来将本文方法与特征矩阵联合近似对角化（joint approximate diagonalization of eigenmatrices，JADE）方法、文献^[6]方法及文献^[11]方法进行对比，采用相似系数作为性能判定指标。在不同的SIR条件下，基于500次蒙特卡罗实验，对比结果如图4所示，其中横轴表示输入的SIR，纵轴表示相似系数。

图4各方法提取目标信号性能对比

Fig.4Target signal extraction performance comparison of different methods

JADE方法是经典的盲信号分离算法，对于本文系统中干扰、回波与直达波信号都有相关性的情况，该方法分离效果不好。文献^[6]方法没有去掉直达波而直接进行回波信号和其他信号的盲分离，即使后续又进行恒模盲均衡算法及时频二维相关算法处理，处理后的强直达波仍会对回波提取性能造成影响，因此回波提取效果也不如本文方法。文献^[11]方法将回波通道中的回波信号投影至由直达波及其时延展开的空间中来消除多径干扰，提取出较为纯净的回波信号，回波信号的提取性能受多径干扰强度的影响较大，因此回波提取效果在SIR较低时不如本文方法。由图4可见，当SIR较低时，本文方法仍具有较好的回波提取性能，仿真结果与理论分析一致。

3.2 实验二：算法目标检测性能分析

受实验条件所限，目前暂时不能获取所需的实测数据，因此利用MATLAB软件仿真实验对本文所提算法和传统方法的目标检测性能进行对比验证，仿真条件如表2所示。

将虚警概率P_fa设定为10^-6，用本文方法恢复出的微弱目标回波信号与直达波信号进行匹配滤波处理，在时延和多普勒平面上提取目标的时延和多普勒峰值信息并评估信噪比，与传统二维匹配滤波方法的检测性能进行对比，结果如图5所示。其中图5（a）和图5（b）分别为传统二维匹配滤波方法和本文方法在时域上的检测性能，图5（c）和图5（d）分别为传统二维匹配滤波方法和本文方法在频域上的检测性能。

表2仿真参数

Tab.2 Simulation parameters

由图5可见，传统二维匹配滤波目标检测方法可以在仿真目标处检测到目标峰值，但目标峰值距平台高度较近且存在虚假峰值；本文方法没有虚假峰值且目标峰值距平台高度相对较远，可以形成理想的目标峰值。利用蒙特卡罗多次试验的方法进行仿真分析，仿真次数为200次，两种方法的目标峰值与平台高度对比如表3所示。

图5本文方法和二维匹配滤波方法目标检测性能对比

Fig.5Target detection performance comparison of the proposed method and the two dimensional matched filtering method

由表3可看出，传统方法中目标峰值与平台高度差值分别为14.91 dB、15.73 dB、14.39 dB，本文方法中目标峰值与平台高度差值分别为20.48 dB、22.64 dB、22.08 dB。该仿真结果表明，传统方法和本文方法都能形成目标峰值，但本文方法性能更优。

表3两种方法目标峰值对比

Tab.3 Comparison of target peaks between the two methods

为了进一步验证本文方法对雷达微弱目标的检测性能，将虚警概率设定为P_fa=10^-6，基于500次蒙特卡罗实验，仿真得到不同的信噪比（signal to noise ratio，SNR）条件下传统二维匹配滤波方法、文献^[6]方法及文献^[11]方法的检测概率P_d的对比曲线，如图6所示。

图6P_fa固定时各种方法P_d随SNR变化曲线对比

Fig.6Comparison of P_d and SNR curves of various methods when P_fa was fixed

由图6可以看出，本文方法在SNR为-15 dB时仍能达到80%的检测概率，在低信噪比下检测目标的性能要优于另外三种方法。将SNR设定为13 dB，基于500次蒙特卡罗实验，仿真得到不同的虚警概率P_fa条件下各种方法的检测概率P_d的ROC曲线，如图7所示，ROC曲线可以从另一个角度来分析这些方法检测目标的性能。

由图7可以看出，当虚警概率P_fa从10^-6变化到1时，本文方法检测目标的性能整体上优于传统二维匹配滤波方法、文献^[6]方法及文献^[11]方法。在虚警概率P_fa较低时，本文方法ROC曲线相对于最优检测曲线来说只有轻微的性能损失，而另外三种方法则有较大的性能损失，无法保持好的检测性能。

图7SNR固定时各种方法P_d随P_fa变化曲线对比

Fig.7Comparison of P_d and P_fa curves of various methods when SNR was fixed

4 结论

对传统GNSS-R外辐射源雷达双通道接收使信号处理过程较为复杂、运算量相对较大，且工程实现上通道结构较为复杂、硬件成本较高的问题，提出一种基于高阶循环累积量的单通道GNSS-R外辐射源雷达微弱目标回波盲检测方法。该方法不需要参考通道，利用PCA的方法及微弱回波信号高阶循环频率与其他信号不同的特性分别将直达波信号及微弱回波信号从单通道接收信号中恢复出来。本文所提方法简化了GNSS-R外辐射源雷达系统的通道结构。仿真结果表明，当SNR较低时本文方法与传统二维匹配滤波目标检测方法及文献方法相比有更好的目标检测性能。

图1GNSS-R外辐射源雷达微弱目标盲检测模型

Fig.1GNSS-R passive radar weak target echo blind detection model

下载: 全尺寸图片

图2单通道微弱目标回波盲检测信号处理流程

Fig.2Signal processing flow chart of single channel weak target echo blind detection

下载: 全尺寸图片

图3源信号及用本文方法恢复出的直达波信号和回波信号

Fig.3Source signals and the direct wave signal and echo signal recovered by the proposed method

下载: 全尺寸图片

图4各方法提取目标信号性能对比

Fig.4Target signal extraction performance comparison of different methods

下载: 全尺寸图片

图5本文方法和二维匹配滤波方法目标检测性能对比

Fig.5Target detection performance comparison of the proposed method and the two dimensional matched filtering method

下载: 全尺寸图片

图6P_fa固定时各种方法P_d随SNR变化曲线对比

Fig.6Comparison of P_d and SNR curves of various methods when P_fa was fixed

下载: 全尺寸图片

图7SNR固定时各种方法P_d随P_fa变化曲线对比

Fig.7Comparison of P_d and P_fa curves of various methods when SNR was fixed

下载: 全尺寸图片

表1相似系数

下载: 全尺寸图片

表2仿真参数

下载: 全尺寸图片

表3两种方法目标峰值对比

下载: 全尺寸图片

图1GNSS-R外辐射源雷达微弱目标盲检测模型

Fig.1GNSS-R passive radar weak target echo blind detection model

图2单通道微弱目标回波盲检测信号处理流程

Fig.2Signal processing flow chart of single channel weak target echo blind detection

图3源信号及用本文方法恢复出的直达波信号和回波信号

Fig.3Source signals and the direct wave signal and echo signal recovered by the proposed method

图4各方法提取目标信号性能对比

Fig.4Target signal extraction performance comparison of different methods

图5本文方法和二维匹配滤波方法目标检测性能对比

Fig.5Target detection performance comparison of the proposed method and the two dimensional matched filtering method

图6P_fa固定时各种方法P_d随SNR变化曲线对比

Fig.6Comparison of P_d and SNR curves of various methods when P_fa was fixed

图7SNR固定时各种方法P_d随P_fa变化曲线对比

Fig.7Comparison of P_d and P_fa curves of various methods when SNR was fixed

表1相似系数

表2仿真参数

表3两种方法目标峰值对比

引用提醒

图(7) / 表(3)

引用本文

温媛媛, 白琳, 尚社, 等. GNSS-R外辐射源雷达微弱目标回波盲检测方法[J]. 国防科技大学学报, 2024, 46(5): 121-130.

复制

WEN Y Y, BAI L, SHANG S, et al. Blind detection method of weak target echo in GNSS-R passive radar[J]. Journal of National University of Defense Technology, 2024, 46(5): 121-130.

Copy

计量

图1GNSS-R外辐射源雷达微弱目标盲检测模型

Fig.1GNSS-R passive radar weak target echo blind detection model

图2单通道微弱目标回波盲检测信号处理流程

Fig.2Signal processing flow chart of single channel weak target echo blind detection

图3源信号及用本文方法恢复出的直达波信号和回波信号

Fig.3Source signals and the direct wave signal and echo signal recovered by the proposed method

图4各方法提取目标信号性能对比

Fig.4Target signal extraction performance comparison of different methods

图5本文方法和二维匹配滤波方法目标检测性能对比

Fig.5Target detection performance comparison of the proposed method and the two dimensional matched filtering method

图6P_fa固定时各种方法P_d随SNR变化曲线对比

Fig.6Comparison of P_d and SNR curves of various methods when P_fa was fixed

图7SNR固定时各种方法P_d随P_fa变化曲线对比

Fig.7Comparison of P_d and P_fa curves of various methods when SNR was fixed

表1相似系数

表2仿真参数

表3两种方法目标峰值对比

WANG T Y, LIU B, QIANG W,et al. GNSS passive radar imaging based on differential aperture[C]//Proceedings of the 18th International Conference on Optical Communications and Networks(ICOCN),2019:1-3.

MA H, ANTONIOU M, PASTINA D,et al. Maritime moving target indication using passive GNSS-based bistatic radar[J]. IEEE Transactions on Aerospace and Electronic Systems,2018,54(1):115-130.

李中余, 黄川, 武俊杰, 等. 基于GNSS的无源雷达海面目标检测技术综述[J]. 雷达科学与技术,2020,18(4):404-416. LI Z Y, HUANG C, WU J J,et al. Overview of maritime target detection techniques using GNSS-based passive radar[J]. Radar Science and Technology,2020,18(4):404-416.(in Chinese)

NASSO I, SANTI F. A centralized approach for ship target detection and localization with multi-transmitters GNSS-based passive radar[C]//Proceedings of International Conference on Radar Systems,2022:202-207.

HE T, CUI L, WANG P B,et al. Non-cooperative air target detection based GNSS bistatic passive radar[C]//Proceedings of 2022 IEEE International Geoscience and Remote Sensing Symposium,2022:2793-2796.

李万春, 田正武, 唐遒, 等. 无参考通道的外辐射源雷达参数估计[J]. 系统工程与电子技术,2016,38(11):2530-2536. LI W C, TIAN Z W, TANG Q,et al. Passive radar parameter estimation without reference channel[J]. Systems Engineering and Electronics,2016,38(11):2530-2536.(in Chinese)

宋海婷. 数字电视外辐射源雷达信号处理技术[D]. 西安: 西安电子科技大学,2018. SONG H T. Signal processing technology of digital TV based passive radar[D]. Xi′an: Xidian University,2018.(in Chinese)

ZHENG H, CHENG S Y, ZHANG M Y. Reference signal reconstruction for GPS-based passive radar based on software radio platform[C]//Proceedings of 2023 IEEE 6th Information Technology, Networking, Electronic and Automation Control Conference(ITNEC),2023:732-736.

CHEN G, TIAN B, GONG J,et al. Reconstruction of passive radar reference signal based on DTMB[C]//Proceedings of the IEEE 2nd International Conference on Information Communication and Signal Processing(ICICSP),2019:170-174.

KRYSIK P, YWEK M. Upper limits of passive radar target detection improvement through removal of noise from reference signal[C]//Proceedings of the Signal Processing Symposium(SPSympo),2019:193-198.

陈刚, 王俊, 王珏, 等. 外辐射源雷达参考信号提纯方法[J]. 系统工程与电子技术,2018,40(1):45-49. CHEN G, WANG J, WANG J,et al. Reference signal purifying method in passive bistatic radar[J]. Systems Engineering and Electronics,2018,40(1):45-49.(in Chinese)

REN X H, LIU N, WANG J H,et al. Direct wave suppression technology of passive radar system[J]. The Journal of Engineering,2019,2019(21):7502-7506.

LIU Y Q, YI J X, WAN X R,et al. Evaluation of clutter suppression in CP-OFDM-based passive radar[J]. IEEE Sensors Journal,2019,19(14):5572-5586.

HE S C, FENG Y, SHAN T. A novel interference suppression method for DTMB-based passive radar[C]//Proceedings of 2021 CIE International Conference on Radar,2021:2799-2803.

ZHANG Q, PENG C, LU Y M,et al. Airborne electromagnetic data levelling using principal component analysis based on flight line difference[J]. Journal of Applied Geophysics,2018,151:290-297.

马愈昭, 陈楠, 熊兴隆. 基于PCA和相位差校正法的风切变预警算法[J]. 系统工程与电子技术,2020,42(1):52-60. MA Y Z, CHEN N, XIONG X L. Wind shear warning algorithm based on PCA and phase difference correction[J]. Systems Engineering and Electronics,2020,42(1):52-60.(in Chinese)

QIN L, LIU S, LONG T,et al. Wind field reconstruction using dimension-reduction of CFD data with experimental validation[J]. Energy,2018,151:272-288.

ŽIVANOVIĆ G D, GARDNER W A. Degrees of cyclostationarity and their application to signal detection and estimation[J]. Signal Processing,1991,22(3):287-297.

HOFMANN J, DELAMOTTE T, KNOPP A. Cyclostationarity-based signal detection in multi-satellite systems[C]//Proceedings of 56th Asilomar Conference on Signals, Systems,and Computers,2022:877-880.

李灯熬, 赵菊敏. 循环平稳理论的盲源分离原理与算法[M]. 北京: 国防工业出版社,2015. LI D A, ZHAO J M. Blind source separation principles and algorithms based on cyclostationary theory[M]. Beijing: National Defense Industry Press,2015.(in Chinese)

首页

期刊介绍

投稿指南

编委会

出版声明

开放获取声明

联系我们

期刊订阅

Rss

AI检索

English

1 单通道GNSS-R外辐射源雷达微弱目标回波盲检测模型

2 单通道GNSS-R外辐射源雷达微弱目标回波盲检测方法

3 仿真分析

4 结论