CN111025257B

CN111025257B - 基于稀疏贝叶斯学习的微动目标高分辨时频图重构方法

Info

Publication number: CN111025257B
Application number: CN201911412455.3A
Authority: CN
Inventors: 白雪茹; 惠叶; 刘思琪; 周峰
Original assignee: Xidian University
Current assignee: Xidian University
Priority date: 2019-12-31
Filing date: 2019-12-31
Publication date: 2023-03-10
Anticipated expiration: 2039-12-31
Also published as: CN111025257A

Abstract

本发明公开了一种基于稀疏贝叶斯学习的微动目标时频图重构方法，主要解决现有技术在低信噪比环境下难以获得微动目标高分辨时频图的问题。其实现方案为：1)获取微动目标的雷达回波信号向量s；2)构建傅里叶字典F；3)从雷达回波信号向量s中获取其每个时间窗内的回波信号s_i，并获取从傅里叶字典F中抽取的每个时间窗对应列构成的字典F_i，利用字典F_i对回波信号s_i进行建模；4)对每个时间窗内的信号s_i及每个时间窗内的时频分布f_i引入稀疏先验；5)基于3)构建的模型及4)中的先验，求解每个时间窗内的时频分布f_i并进行拼接，得到微动目标的时频图重构结果。本发明能在低信噪比下获得弹道目标、微动目标聚焦良好的时频图，可用于复杂环境微动特征提取与识别。

Description

基于稀疏贝叶斯学习的微动目标高分辨时频图重构方法

技术领域

本发明属于雷达信号处理技术领域，特别涉及一种微动目标高分辨时频图重构方法，可用于复杂观测条件下的微动目标高分辨成像、特征提取与识别。

背景技术

微动是指目标或目标上的某些结构存在的独立于主体运动的振动或旋转，例如直升机叶片的旋转、弹道目标的自旋、进动和章动，以及人体四肢的摆动等。微动会在目标主体运动对应的主多普勒周围产生边带，即产生微多普勒效应。由于微动目标参数化特征提取方法通常需要建立复杂的参数化模型，并且求解运算量很大。因此基于联合时频JTF分析的非参数化方法被广泛用于微动目标成像与特征提取，该类方法对微动形式具有鲁棒性，计算效率较高。然而，在低信噪比的观测环境下，利用传统时频分析方法所得时频图将出现虚假点，从而影响准确的目标特征提取。

针对低信噪比等观测环境，可采用稀疏信号重构技术实现有效的JTF分析。例如：白雪茹，周峰在其发表的论文“Radar imaging of micromotion targets from corrupteddata”(IEEE Trans.Aerosp.&Electron.Syst..vol.52,no.6,pp.2789-2802,Dec.2016.)中，公开了一种在低信噪比条件下获得微动目标高质量时频图的方法。该方法将JTF重构问题转化为L1范数最优化问题，从而获得微动目标高质量时频图，进而实现目标高分辨成像。但是该类算法难以确定正则化参数。

发明内容

本发明目的在于提出一种基于稀疏贝叶斯学习的微动目标高分辨时频图重构方法，以解决在上述现有技术在低信噪比的观测环境下难以获得微动目标高质量时频图的问题。

本发明的技术方案是：先利用非参数化字典对微动目标回波信号进行JTF域建模，然后构建基于伽马-复高斯先验的概率模型，进而采用变分推断求解模型中的未知系数，得到高分辨时频图,具体实现步骤包括如下：

(1)发射线性调频信号进行解线频调处理，并接收低信噪比环境下微动目标维数为N×1的雷达回波信号向量s，其中，N为采样点数；

(2)以e^j2πmn/N为元素，构造维数为N×N的傅里叶字典F，其中，e⁽·⁾表示以自然常数为底的指数操作，j表示虚数单位，m表示傅里叶字典F的行序号，n表示傅里叶字典F的列序号，行序号m与列序号n的取值范围均为[1:N]；

(3)获取每个时间矩形滑窗内的信号s_i：

设L为每个个滑窗的窗长，设s_i表示第i个滑窗内的回波信号，其维度为L×1，i∈[1,N]；

对每个时间矩形滑窗内的信号s_i进行建模：

s_i＝G_if_i+n_i

其中，G_i为矩阵F_iH_i的伪逆矩阵，其维度为L×N；H_i表示第i个滑窗窗函数对应的对角矩阵，其维度为L×L，F_i表示从矩阵F中抽取第i个滑窗内对应的列所构成的第i个滑窗对应的傅里叶字典，其维度为N×L；n_i为噪声向量，维度为L×1，f_i为第i个滑窗内的时频分布，f_i＝[f_i,1,...,f_i,n,...,f_i,N]^T,其中f_i,n为f_i的第n个元素，n＝1,...,N,f_i维度为N×1；；

(4)假设每个时间矩形滑窗内的信号s_i及每个时间矩形滑窗内的时频分布f_i服从伽马-复高斯先验，则s_i的先验概率密度函数p(s_i|α_i,f_i)及f_i的先验概率密度函数p(f_i|λ_i)的表达式分别为：

p(f_i|λ_i)＝CN(f_i|0,Σ_i)

其中，CN(·)表示复高斯分布，(·)^H表示矩阵的共轭转置，α_i是第i个滑窗内的噪声精度参数，每个滑窗内对应的α_i独立，

α_i服从伽马分布p(α_i)，λ_i是第i个滑窗内时频分布f_i的分布参数，λ_i服从伽马分布p(λ_i)，p(α_i)及p(λ_i)的表达式分别为：

其中，Gam(·)表示伽马分布，c_i，b_i，u_i,n和v_i,n为伽马分布的不同参数；

(5)依据(3)构建的模型及(4)中的先验概率密度函数，利用变分推断方法求解每个时间矩形滑窗内的时频分布f_i：

(5a)将阈值T的值设置为10-⁶，最大迭代次数M的值设置为500；

(5b)令滑窗的初始序号为i＝1；

(5c)初始化参数：

令初始迭代次数为k＝1，将参数c_i，b_i的初始值设置为

将参数u_i,n的初始值设置为

参数v_i,n的初始值设置为

第i个滑窗内的噪声精度参数α_i初始值

及第i个滑窗内时频分布f_i的分布参数λ_i的初始值

根据回波信噪比的不同而设置；

(5d)计算第k次迭代过程中第i个滑窗内时频分布f_i的近似后验分布q^(k)(f_i)：

(5e)计算与第i个滑窗内的时频分布f_i的三个相关期望：

和

(5f)分别计算第i个滑窗内噪声精度参数α_i的第k次迭代的值

和第i个滑窗内时频分布f_i的分布参数λ_i的第k次迭代的值

(5g)判断当前迭代次数k＞1是否成立，若是，执行(5h)；否则，令k＝k+1，返回(5d)；

(5h)判断

或k＞M是否成立，若是，则第i个滑窗内的时频分布f_i为

并执行步骤(5i)；否则，令k＝k+1，返回(5d)，其中，||·||₂表示2-范数操作；

(5i)判断i＝N是否成立，若是，得到每个时间矩形滑窗内的时频分布f_i；否则，令i＝i+1，并执行步骤(5c)；

(6)将每个时间矩形滑窗内的时频分布f_i，按行拼接，得到一个N×N的矩阵，即为微动目标的时频图重构结果。

本发明由于构建了基于伽马-复高斯先验的概率模型，采用变分推断求解模型中的稀疏系数，具有以下优点：

1)在低信噪比情况下，能够得到稀疏且高分辨的时频图，克服了现有技术在低信噪比条件下重构时频图存在虚假点的不足。

2)避免了现有基于L1范数优化方法在JTF重构时需要预先确定正则化参数的不足，简化了高分辨时频图的重构过程。

附图说明

图1为本发明的实现流程图；

图2为本发明的仿真结果图。

具体实施方式

下面结合附图对本发明的实施例和效果作进一步详细描述。

参照图1，本实例基于贝叶斯学习的微动目标高分辨时频图重构方法，包括如下步骤：

步骤1，获取微动目标的雷达回波信号向量s。

发射线性调频信号进行解线频调处理，并接收低信噪比环境下微动目标的回波信号，对回波信号进行采样，取采样点数为N，得到维数为N×1的雷达回波信号向量s。

步骤2，构建傅里叶字典F。

以e^j2πmn/N为元素，构造维数为N×N的傅里叶字典F：

其中，e^(·)表示以自然常数为底的指数操作，j表示虚数单位，m表示傅里叶字典F的行序号，n表示傅里叶字典F的列序号，行序号m与列序号n的取值范围均为[1:N]。

步骤3，获取每个时间矩形滑窗内的信号s_i，并对其进行建模。

设L为各个滑窗的窗长，其中滑窗可为矩形窗，汉明窗，三角窗，本实例采用的是矩形窗，窗长设置为L_i＝31；

设s_i表示第i个滑窗内的回波信号，其维度为L×1，i∈[1,N]；

对每个时间矩形滑窗内的信号s_i进行建模：

s_i＝G_if_i+n_i

其中，G_i为矩阵F_iH_i的伪逆矩阵，其维度为L×N；H_i表示第i个滑窗窗函数对应的对角矩阵，其维度为L×L，F_i表示从矩阵F中抽取第i个滑窗内对应的列所构成的第i个滑窗对应的傅里叶字典，其维度为N×L；n_i为噪声向量，维度为L×1，f_i为第i个滑窗内的时频分布，f_i＝[f_i,1,...,f_i,n,...,f_i,N]^T,其中f_i,n为f_i的第n个元素，n＝1,...,N，,f_i的维度为N×1。

步骤4，求解每个时间矩形滑窗内的时频分布f_i。

4.1)假设每个时间矩形滑窗内的信号s_i及每个时间矩形滑窗内的时频分布f_i服从伽马-复高斯先验，则s_i的先验概率密度函数p(s_i|α_i,f_i)及f_i的先验概率密度函数p(f_i|λ_i)的表达式分别为：

p(f_i|λ_i)＝CN(f_i|0,Σ_i)

4.2)依据上述先验概率密度函数和(3)构建的模型，利用变分推断方法求解每个时间矩形滑窗内的时频分布f_i：

(4.2.1)将阈值T的值设置为10-⁶，最大迭代次数M的值设置为500；

(4.2.2)令滑窗的初始序号为i＝1；

(4.2.3)初始化参数：

令初始迭代次数为k＝1，将参数c_i，b_i的初始值设置为

将参数u_i,n的初始值设置为

参数v_i,n的初始值设置为

根据回波信噪比的不同设置第i个滑窗内的噪声精度参数α_i初始值

及第i个滑窗内时频分布f_i的分布参数λ_i的初始值

即在回波信噪比为5dB的情况下α_i的初始值设置为

λ_i的初始值

中的第d个元素

的初始值设置为

(4.2.4)计算第k次迭代过程中第i个滑窗内时频分布f_i的近似后验分布q^(k)(f_i)：

其中，

为近似后验分布q^(k)(f_i)的协方差矩阵，

为近似后验分布q^(k)(f_i)的均值；

(4.2.5)计算与第i个滑窗内的时频分布f_i的三个相关期望：

和

其中，

为近似后验分布q^(k)(f_i)的均值，

为

的第n个元素，

为近似后验分布q^(k)(f_i)的协方差矩阵，

为矩阵

的第n行第n列元素；

(4.2.6)分别计算第i个滑窗内噪声精度参数α_i的第k次迭代的值

和第i个滑窗内时频分布f_i的分布参数λ_i的第k次迭代的值

其中，

(4.2.7)判断当前迭代次数k＞1是否成立，若是，则执行(4.2.8)；否则，令k＝k+1，返回(4.2.4)；

(4.2.8)判断

或k＞M是否成立，若是，则第i个滑窗内的时频分布f_i为

并执行步骤(4.2.9)；否则，令k＝k+1，返回(4.2.4)，其中，||·||₂表示2-范数操作；

(4.2.9)判断i＝N是否成立，若是，得到每个时间矩形滑窗内的时频分布f_i；否则，令i＝i+1，返回步骤(4.2.3)。

步骤6，将得到的各滑窗内的时频分布f_i，i＝1,...,N按行拼接，得到一个N×N的矩阵，该矩阵为微动目标的时频图重构结果。

以下结合仿真实验，对本实例的技术效果作进一步说明。

1.仿真条件：

本实例的仿真实验采用工作在C波段的雷达，对应载频为10GHZ，脉冲重复频率为500Hz，采样点数取N＝256。

目标由4个散射点组成，目标自旋频率为0.5Hz。

2.仿真内容：

仿真1，对信噪比为5dB的回波信号做STFT，将结果绘制为目标的时频仿真图，如图2(a)所示，其中横坐标表示慢时间，纵坐标表示目标的多普勒频率。由图2(a)可以看出，用现有方法在低信噪比条件下，所得的时频图会出现虚假点，影响准确的目标特征提取。

仿真2，利用本发明对信噪比为5dB的回波信号进行时频图重构，绘制仿真结果图，如图2(b)所示，横坐标表示慢时间，纵坐标表示目标的多普勒频率。图2(b)可以看出，利用本发明可以得到高分辨时频图。

由上述仿真结果表明，采用本发明可以在低信噪比情况下获得高分辨时频图。与现有技术相比，本发明避免了现有基于L1范数优化方法在JTF重构时需要预先确定正则化参数的不足，解决了在低信噪比的观测环境下难以获得微动目标高分辨时频图的问题。

Claims

1.一种基于稀疏贝叶斯学习的微动目标高分辨时频图重构方法，其特征在于，包括如下步骤：

(2)以e^j2πmn/N为元素，构造维数为N×N的傅里叶字典F，其中，e^(·)表示以自然常数为底的指数操作，j表示虚数单位，m表示傅里叶字典F的行序号，n表示傅里叶字典F的列序号，行序号m与列序号n的取值范围均为[1:N]；

(3)获取每个时间矩形滑窗内的信号s_i：

设L为每个滑窗的窗长，设s_i表示第i个滑窗内的回波信号，其维度为L×1，i∈[1,N]；

对每个时间矩形滑窗内的信号s_i进行建模：

s_i＝G_if_i+n_i

其中，G_i为矩阵F_iH_i的伪逆矩阵，其维度为L×N；H_i表示第i个滑窗窗函数对应的对角矩阵，其维度为L×L，F_i表示从矩阵F中抽取第i个滑窗内对应的列所构成的第i个滑窗对应的傅里叶字典，其维度为N×L；n_i为噪声向量，维度为L×1，f_i为第i个滑窗内的时频分布，f_i＝[f_i,1,...,f_i,n,...,f_i,N]^T,其中f_i,n为f_i的第n个元素，n＝1,...,N,f_i维度为N×1；