基于SET-STIAA的空间微动目标时频分析方法

摘要

本发明提供的一种基于SET‑STIAA的空间微动目标时频分析方法，将同步提取变换与短时迭代自适应算法相结合，利用IAA算法，对每个时刻的滑窗内的两个回波信号计算频率成分，获得每个时刻的滑窗内每个回波信号的第一频率曲线函数；使用同步提取变换方法在第一频率曲线函数中选择瞬时相位对应的回波信号频点作为该回波信号估计的瞬时频率；遍历每一个瞬时频率点，将该瞬时频率点的频率与该频点处的真实频率进行比较，以剔除频谱泄露的瞬时频率点从而确定回波信号的频率。本发明仅提取出瞬时相位处的时频点，达到提高时频分辨精度的效果。并且只需要更少的滑窗点数，有效缓解数据长度较短且微多普勒曲线距离较近时不能分辨出目标的缺点。

权利要求

1.一种基于SET‑STIAA的空间微动目标时频分析方法，其特征在于，包括：将雷达系统中的两个回波信号进行滑窗处理；

利用IAA算法，对每个时刻的滑窗内的两个回波信号计算频率成分，获得每个时刻的滑窗内每个回波信号的第一频率曲线函数；

针对每个回波信号的第一频率曲线函数，使用同步提取变换方法，在第一频率曲线函数中选择瞬时相位对应的回波信号频点作为该回波信号估计的瞬时频率；

将每个回波信号估计的瞬时频率组成时频曲线图；

在每个回波信号的时频曲线图中遍历每一个瞬时频率点，将该瞬时频率点的频率与该频点处的真实频率进行比较，以剔除频谱泄露的瞬时频率点，获得每个回波信号的频率曲线图以确定回波信号的频率；

其中，两个回波信号同一频点的频率值之差大于2倍的频率间隔；

所述使用同步提取变换方法在所述第一频率曲线函数中选择瞬时相位对应的回波信号频点作为回波信号估计的瞬时频率包括：使用同步提取变换方法对所述IAA算法进行变换，获得第二频率曲线函数；

使用回波信号的泰勒展开式以及第一频率曲线函数，对所述第二频率曲线函数进行转化，获得求取所述瞬时频率的表达式；

使用所述瞬时频率的表达式在所述第一频率曲线函数中选择瞬时相位对应的回波信号频点作为回波信号估计的瞬时频率。

2.根据权利要求1所述的目标时频分析方法，其特征在于，在所述利用IAA算法，对每个滑窗内的每个回波信号计算频率成分，获得在滑窗内的每个回波信号的第一频率曲线函数之前，所述目标时频分析方法还包括：将所述回波信号进行离散变换，获得离散变换后的回波信号；

将所述离散变换后的回波信号进行泰勒展开，获得回波信号的泰勒展开式。

3.根据权利要求2所述的目标时频分析方法，其特征在于，所述回波信号表示为：离散变换后的回波信号表示为：

回波信号的泰勒展开式为：

其中，y(t)表示回波信号，t表示时间，A(t)表示信号的幅值函数， 表示信号的相位函数，y(n)表示离散变换后的回波信号，A(n)表示离散变换后的回波信号的幅值， 表示离散变换后回波信号的相位，n表示离散时刻，y(n1)表示回波信号在n1处进行泰勒展开，表示回波信号的真实频率，i表示虚数单位。

4.根据权利要求3所述的目标时频分析方法，其特征在于，所述第一频率曲线函数表示为：所述第二频率曲线表示为：

瞬时频率的表达式为：

H

其中，G(K,w)代表K时刻滑窗在w频点出的幅值，a(w)代表回波信号的傅里叶积，a (w)代‑1表回波信号逆傅里叶积的共轭转置，RK 代表回波信号相关函数的逆矩阵，yk(n)代表在k时刻的回波信号，yk(n)＝y(n)g(n‑k)，g(n‑k)表示将窗函数延迟n时刻，g(k)表示窗函数，P(K,w)表示k+1时刻滑窗在ω频点出的幅值，a(w)代表信号的傅里叶积分，yk+1(n)表示后一个脉冲时刻的回波信号，K表示第K个时刻的滑窗，ω(k,w)表示估计的瞬时频率，T表示观测的总时间长度，N表示观测时间内的脉冲数。

5.根据权利要求4所述的目标时频分析方法，其特征在于，所述在每个回波信号的时频曲线图中遍历每一个瞬时频率点，将该瞬时频率点的频率与该频点处的真实频率进行比较，以剔除频谱泄露的瞬时频率点，获得每个回波信号的频率曲线图包括：使用SET算法，在每个回波信号的时频曲线图中遍历每一个瞬时频率点，将该瞬时频率点的频率与该频点处的真实频率进行比较，以剔除频谱泄露的瞬时频率点，获得每个回波信号的频率曲线图；

其中，SET算法的表达式为：Ts(k,γ)＝δ(γ‑ω(k,w))G(k,w)，δ(γ‑ω(k,w))表示选择函数，γ表示频点处的真实频率值。

说明书

基于SET‑STIAA的空间微动目标时频分析方法

技术领域

[0001] 本发明属于雷达技术领域，具体涉及一种基于SET‑STIAA的空间微动目标时频分析方法。

背景技术

[0002] 空间微动目标的微多普勒特征提取是空间椎体目标识别中最为关键的步骤之一，而特征提取的重心是对目标的参数来进行估计，所以时频分步逐渐成为研究人员提取微动多普勒变换曲线最为直观的反应。

[0003] 现有的时频分析方法主要有以下几种方法：

[0004] (1)基于同步提取变换的短时傅里叶变换(SET‑STFT)，是最常用的也是最为简单的时频分析方法，它是由傅里叶变换发展而来。SET‑STFT在传统傅里叶变换的基础上采用时间滑窗将一系列非平稳信号分解为若干短时有限段，在每一小段内认为该信号是平稳的，对滑窗截取的每一小段原始信号进行加窗处理，以降低因信号截取造成的频谱泄露效应，然后对加窗信号做傅里叶变换。但是基于同步提取变换的短时傅里叶变换不适用于距离过近的两条曲线，并且在曲线的交叉处分辨不开，同时SET‑STFT算法所需要的短时滑窗点数相对较高，在计算一段连续的信号时需要舍弃一个滑窗的点数，不适用于观测时间较短的信号。

[0005] (2)短时迭代自适应算法(STIAA)，是将非参数的迭代自适应算法应用于时频分析领域，发展了适用于波达方向估计领域的短时迭代自适应算法，并且均达到了时频高分辨的效果。STIAA算法是在IAA算法的基础上发展而来的，对一定长度的原始离散信号做类似于STFT算法的滑窗处理，对每一个滑窗内的信号应用IAA算法。但是STIAA还是与理想的时频分辨率有一定的差距，不能满足要求。

发明内容

[0006] 为了解决现有技术中存在的上述问题，本发明提供了一种基于SET‑STIAA的空间微动目标时频分析方法。本发明要解决的技术问题通过以下技术方案实现：

[0007] 本发明提供的一种基于SET‑STIAA的空间微动目标时频分析方法包括：

[0008] 将雷达系统中的两个回波信号进行滑窗处理；

[0009] 利用IAA算法，对每个时刻的滑窗内的两个回波信号计算频率成分，获得每个时刻的滑窗内每个回波信号的第一频率曲线函数；

[0010] 针对每个回波信号的第一频率曲线函数，使用同步提取变换方法，在第一频率曲线函数中选择瞬时相位对应的回波信号频点作为该回波信号估计的瞬时频率；

[0011] 将每个回波信号估计的瞬时频率组成时频曲线图；

[0012] 在每个回波信号的时频曲线图中遍历每一个瞬时频率点，将该瞬时频率点的频率与该频点处的真实频率进行比较，以剔除频谱泄露的瞬时频率点，获得每个回波信号的频率曲线图以确定回波信号的频率；

[0013] 其中，两个回波信号同一频点的频率值之差大于2倍的频率间隔。

[0014] 可选的，在所述利用IAA算法，对每个滑窗内的每个回波信号计算频率成分，获得在滑窗内的每个回波信号的第一频率曲线函数之前，所述目标时频分析方法还包括：

[0015] 将所述回波信号进行离散变换，获得离散变换后的回波信号；

[0016] 将所述离散变换后的回波信号进行泰勒展开，获得回波信号的泰勒展开式。

[0017] 其中，所述回波信号表示为：

[0018]

[0019] 离散变换后的回波信号表示为：

[0020]

[0021] 回波信号的泰勒展开式为：

[0022]

[0023] 其中，y(t)表示回波信号，t表示时间，A(t)表示信号的幅值函数， 表示信号的相位函数，y(n)表示离散变换后的回波信号，A(n)表示离散变换后的回波信号的幅值，表示离散变换后回波信号的相位，n表示离散时刻，n1表示回波信号在n1处进行泰勒展开，表示回波信号的真实频率，i表示虚数单位。

[0024] 可选的，所述使用同步提取变换方法在所述第一频率曲线函数中选择瞬时相位对应的回波信号频点作为回波信号估计的瞬时频率包括：

[0025] 使用同步提取变换方法对所述IAA算法进行变换，获得第二频率曲线函数；

[0026] 使用回波信号的泰勒展开式以及第一频率曲线函数，对所述第二频率曲线函数进行转化，获得求取所述瞬时频率的表达式；

[0027] 使用所述瞬时频率的表达式在所述第一频率曲线函数中选择瞬时相位对应的回波信号频点作为回波信号估计的瞬时频率。

[0028] 其中，所述第一频率曲线函数表示为：

[0029]

[0030] 所述第二频率曲线表示为：

[0031] 瞬时频率的表达式为：

[0032] 其中，G(K,w)代表K时刻滑窗在w频点出的幅值，a(w)代表回波信号的傅里叶积，aH‑1(w)代表回波信号逆傅里叶积的共轭转置，RK 代表回波信号相关函数的逆矩阵，yk(n)代表在k时刻的回波信号，yk(n)＝y(n)g(n—k)，g(n—k)表示将窗函数延迟n时刻，g(k)表示窗函数， P(K,w)表示k+1时刻滑窗在ω频点出的幅值，a
(w)代表信号的傅里叶积分，yk+1(n)表示后一个脉冲时刻的回波信号，K表示第K个时刻的滑窗，ω(k,w)表示估计的瞬时频率，T表示观测的总时间长度，N表示观测时间内的脉冲数。

[0033] 可选的，所述在每个回波信号的时频曲线图中遍历每一个瞬时频率点，将该瞬时频率点的频率与该频点处的真实频率进行比较，以剔除频谱泄露的瞬时频率点，获得每个回波信号的频率曲线图包括：

[0034] 使用SET算法，在每个回波信号的时频曲线图中遍历每一个瞬时频率点，将该瞬时频率点的频率与该频点处的真实频率进行比较，以剔除频谱泄露的瞬时频率点，获得每个回波信号的频率曲线图；

[0035] 其中，SET算法的表达式为：Ts(k,γ)＝δ(γ—ω(k,w))G(k,w)，δ(γ—ω(k,w))表示选择函数，γ表示频点处的真实频率值。

[0036] 1、本发明提供的一种基于SET‑STIAA的空间微动目标时频分析方法，将同步提取变换与短时迭代自适应算法相结合，使用滑窗获取雷达系统中的两个回波信号；利用IAA算法，对每个时刻的滑窗内的两个回波信号计算频率成分，获得每个时刻的滑窗内每个回波信号的第一频率曲线函数；使用同步提取变换方法，第一频率曲线函数中选择瞬时相位对应的回波信号频点作为该回波信号估计的瞬时频率；遍历每一个瞬时频率点，将该瞬时频率点的频率与该频点处的真实频率进行比较，以剔除频谱泄露的瞬时频率点，从而确定回波信号的频率，本发明仅提取出瞬时相位处的时频点，达到提高时频分辨精度的效果。

[0037] 2、本发明在对微动目标的时频分析的过程中，只需要更少的滑窗点数，有效缓解了数据长度较短且微多普勒曲线距离较近时不能分辨出目标的缺点。

[0038] 以下将结合附图及实施例对本发明做进一步详细说明。

附图说明

[0039] 图1为本发明提供的一种基于SET‑STIAA的空间微动目标时频分析方法的流程图；

[0040] 图2为本发明对无噪声电磁仿真数据进行分析的结果图；

[0041] 图3为本发明对15dB噪声电磁仿真数据进行分析的结果图。

具体实施方式

[0042] 下面结合具体实施例对本发明做进一步详细的描述，但本发明的实施方式不限于此。

[0043] 实施例一

[0044] 如图1所示，本发明提供的一种基于SET‑STIAA的空间微动目标时频分析方法包括：

[0045] S1，将雷达系统中的两个回波信号进行滑窗处理；

[0046] 可以理解，在进行滑窗处理前，需要从雷达系统中获取两个回波信号，获取方式可以根据现有的获取方式，此处不再赘述。

[0047] S2，利用IAA算法，对每个时刻的滑窗内的两个回波信号计算频率成分，获得每个时刻的滑窗内每个回波信号的第一频率曲线函数；

[0048] 其中，两个回波信号同一频点的真实频率之差大于2倍的频率间隔。

[0049] 对于两个信号模型，只有满足以下条件才能应用SET算法达到提高时频分辨率的目的：

[0050]

[0051] 其中k∈[2,3,...,K]代表信号分量的总个数，Δw表示多个信号之间的频率间隔，表示信号的真实频率。

[0052] S3，针对每个回波信号的第一频率曲线函数，使用同步提取变换方法，第一频率曲线函数中选择瞬时相位对应的回波信号频点作为该回波信号估计的瞬时频率；

[0053] S4，将每个回波信号估计的瞬时频率组成时频曲线图；

[0054] S5，在每个回波信号的时频曲线图中遍历每一个瞬时频率点，将该瞬时频率点的频率与该频点处的真实频率进行比较，以剔除频谱泄露的瞬时频率点，获得每个回波信号的频率曲线图；

[0055] 本发明提供的一种基于SET‑STIAA的空间微动目标时频分析方法，将同步提取变换与短时迭代自适应算法相结合，利用IAA算法，对每个时刻的滑窗内的两个回波信号计算频率成分，获得每个时刻的滑窗内每个回波信号的第一频率曲线函数；使用同步提取变换方法在第一频率曲线函数中选择瞬时相位对应的回波信号频点作为该回波信号估计的瞬时频率；遍历每一个瞬时频率点，将该瞬时频率点的频率与该频点处的真实频率进行比较，以剔除频谱泄露的瞬时频率点从而确定回波信号的频率。本发明仅提取出瞬时相位处的时频点，达到提高时频分辨精度的效果。并且只需要更少的滑窗点数，有效缓解数据长度较短且微多普勒曲线距离较近时不能分辨出目标的缺点。

[0056] 作为本发明一种可选的实施方式，在利用IAA算法，对每个滑窗内的每个回波信号计算频率成分，获得在滑窗内的每个回波信号的第一频率曲线函数之前，目标时频分析方法还包括：

[0057] 步骤一：将回波信号进行离散变换，获得离散变换后的回波信号；

[0058] 步骤二：将离散变换后的回波信号进行泰勒展开，获得回波信号的泰勒展开式。

[0059] 假设回波信号y(t)的模型满足下式：

[0060]

[0061] 同时假定存在ε足够小，对于任意的t使得|A'(t)＜ε|并且满足 其中A(t)表示信号的幅值函数， 表示信号的相位函数。

[0062] 离散变换后的回波信号表示为：

[0063]

[0064] 回波信号的泰勒展开式为：

[0065]

[0066] 其中，y(t)表示回波信号，t表示时间，A(t)表示信号的幅值函数， 表示信号的相位函数，y(n)表示离散变换后的回波信号，A(n)表示离散变换后的回波信号的幅值，表示离散变换后回波信号的相位，n表示离散时刻，n1表示回波信号在n1处进行泰勒展开，表示回波信号的真实频率，i表示虚数单位。

[0067] 作为本发明一种可选的实施方式，使用同步提取变换方法在第一频率曲线函数中选择瞬时相位对应的回波信号频点作为回波信号估计的瞬时频率包括：

[0068] 步骤一：使用同步提取变换方法对IAA算法进行变换，获得第二频率曲线函数；

[0069] 短时迭代自适应算法STIAA公式定义如下：

[0070]

[0071] 所以第一频率曲线函数表示为：

[0072]

[0073] 其中，G(K,w)代表K时刻滑窗在w频点出的幅值，a(w)代表信号的傅里叶积，aH(w)‑1代表信号逆傅里叶积的共轭转置，RK 代表输入信号相关函数的逆矩阵，yk(n)代表在k时刻的输入信号。

[0074] 步骤二：使用回波信号的泰勒展开式以及第一频率曲线函数，对第二频率曲线函数进行转化，获得求取瞬时频率的表达式；

[0075] 步骤三：使用瞬时频率的表达式在第一频率曲线函数中选择瞬时相位对应的回波信号频点作为回波信号估计的瞬时频率。

[0076] 第二频率曲线表示为：

[0077] 可以推导出估计瞬时频率的公式，对于任意的G(k,w)≠0并且P(k,w)≠0，瞬时频率估计ω(k,w)为：

[0078] 其中，G(K,w)代表K时刻滑窗在w频点出的幅值，a(w)代表回波信号的傅里叶积，aH‑1(w)代表回波信号逆傅里叶积的共轭转置，RK 代表回波信号相关函数的逆矩阵，yk(n)代表在k时刻的回波信号，yk(n)＝y(n)g(n‑k)，g(n‑k)表示将窗函数延迟n时刻，g(k)表示窗函数， P(K,w)表示k+1时刻滑窗在ω频点出的幅值，a(w)
代表信号的傅里叶积，yk+1(n)表示后一个脉冲时刻的回波信号，K表示第K个时刻的滑窗，ω(k,w)表示估计的瞬时频率，T表示观测的总时间长度，N表示观测时间内的脉冲数(慢时间)。

[0079] 作为本发明一种可选的实施方式，在每个回波信号的时频曲线图中遍历每一个瞬时频率点，将该瞬时频率点的频率与该频点处的真实频率进行比较，以剔除频谱泄露的瞬时频率点，获得每个回波信号的频率曲线图包括：

[0080] 使用SET算法，在每个回波信号的时频曲线图中遍历每一个瞬时频率点，将该瞬时频率点的频率与该频点处的真实频率进行比较，以剔除频谱泄露的瞬时频率点，获得每个回波信号的频率曲线图；

[0081] 其中，SET算法的表达式为：Ts(k,γ)＝δ(γ‑ω(k,w))G(k,w)，δ(γ‑ω(k,w))表示选择函数，γ表示频点处的真实频率值。

[0082] 由于一个时频曲线做滑窗处理之后会有频谱泄露的问题，时频图会变宽，导致时频曲线分辨率会降低，为了去掉频谱泄露的部分，留下真实的频谱成分。

[0083] 证明，令 并且ε足够小的情况下，对于任意的 并且 可以得到如下近似；

[0084]

[0085] 通过上式可以看出来瞬时频率的估计值ω(k,w)可以非常接近真实的信号瞬时频率。因此将用于频率提取以提高时频分辨精度的SET算法的定义如下：

[0086] Ts(k,γ)＝δ(γ‑ω(k,w))G(k,w)，

[0087] 表示信号的真实频率，ε表示任意小的数。

[0088] 其中，选择函数对计算出来的频率进行选择，从1开始遍历到最后一个频点。

[0089] 两个信号分量的模型进行如下分析：

[0090] y(n)＝y1(n)+y2(n)

[0091] 瞬时频率估计为：

[0092]

[0093] 瞬时频率估计ω(k,w)是一个复数，而真实的瞬时频率是一个实数，对瞬时频率估计选择取其实部。如果信号y1和信号y2满足两个回波信号同一频点的真实频率之差大于2倍的频率间隔，则可以推导出下式：

[0094]

[0095] 表示信号y1的真实频率， 表示信号y2的真实频率， 表示在信号y1的真实频率所在的区间， 表示在信号y2真实频率所在的区间。

[0096] 将上式代入式Ts(k,γ)＝δ(γ‑ω(k,w))G(k,w)中，即可得到最终结果‑基于同步压缩变换的高分辨时频分析结果。

[0097] 下面结合仿真实验对本发明的效果做进一步描述。

[0098] 1.仿真实验条件。

[0099] 本发明的仿真实验的硬件平台为：处理器为Inter(R)Core(TM)i7‑4790CPU主频为3.60GHZ，内存32GB。

[0100] 本发明的仿真实验的软件平台为：CST STUDIO SUITE 2019和MATLAB R2018b。

[0101] 2.仿真内容。

[0102] 本发明的仿真实验是采用本发明和两个现有技术(短时傅里叶变换STFT时频分析方法、短时迭代自适应时频分析方法)在分别对如图2(a)所示的空间微动目标模型的回波信号进行时频分析，获得三种技术对应的时频图。

[0103] 用于仿真实验的空中微动目标模型的主要参数设置如表1所示。

[0104] 表1空中微动目标模型主要参数表

[0105]参数 数值
底部半径/m 0.44
顶部半径/m 0.0488
目标高度/m 1.77
自旋角频率/HZ 0
进动角频率/HZ 1
章动角频率/HZ 0

[0106] 仿真数据参数设置如下：

[0107] 滑窗长度16点，超分辨倍数20，进动频率1Hz，脉冲重复频率800Hz。

[0108] 图2为无噪声电磁仿真结果，图3为15dB噪声电磁仿真结果。

[0109] 3.仿真结果分析

[0110] 下面结合图2、图3的仿真图对本发明的效果作进一步的描述

[0111] 图2中的(a)为短时傅里叶变换方法时频分析结果，时频分辨精度很差，对于较近的两条交叉曲线无法分辨开来。图2中的(b)为短时迭代自适应方法时频分析结果，虽然对时频分辨率精度有所提高，能够分辨出距离较近的两条交叉曲线，但是能看见明显的毛刺，曲线刻画依旧不够精细。图2中的(c)为基于同步提取变换的短时迭代自适应方法的时频分析结果，能够明显分辨出距离较近的两条交叉曲线，与图2中的(b)相比，时频分辨精度明显提高，曲线刻画更为精细，更接近理想的时频曲线。

[0112] 图3是对带有15dB噪声的电磁仿真数据进行时频分析的结果，图3中的(a)、(b)、(c)分别为短时傅里叶变换方法时频分析结果、短时迭代自适应方法时频分析结果、基于同步提取变换的短时迭代自适应方法的时频分析结果，可以清楚的看出基于同步提取变换的短时迭代自适应方法的时频分析结果分辨精度更高、曲线刻画更精细。可见基于同步提取变换的短时迭代自适应的时频分析方法抗噪性能较好。

[0113] 以上仿真实验表明：本发明的时频分析方法，能够得到时频分辨率较高的时频图，对于距离较近的曲线，可以很清晰的将两条交叉曲线分辨开，并且抗噪声性能相对较好，分辨精度得到了极大的提升。以上内容是结合具体的优选实施方式对本发明所作的进一步详细说明，不能认定本发明的具体实施只局限于这些说明。对于本发明所属技术领域的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干简单推演或替换，都应当视为属于本发明的保护范围。