CN115795233A - 频域数据的延时项的提取方法、装置及存储介质 - Google Patents

Abstract

本发明涉及一种频域数据的延时项的提取方法、装置及存储介质。本发明将传统的Gabor变换的时‑频域分解方法中的频率积分上下限近似为无穷，再运用高斯分布的概率密度积分公式对计算式进行化简，直接使用输入的频域数据得到具有瞬时功率特性的时域结果，大大提高了运算效率。

Description

频域数据的延时项的提取方法、装置及存储介质

技术领域

本发明属于射频信号处理领域，从时频分析领域发展而来，更具体的说涉及的是对含有延时的频域数据进行延时提取的快速方法。

背景技术

在射频多端口设备的仿真中，通常利用频域数据，例如用S参数来代表无源网络。S参数是散射参数的缩写，其描述了传输通道的频域特性，例如：信号的反射、串扰、损耗、延时等信息，都可以通过S参数找到相关的信息，S参数在射频、微波的应用都很广泛，它是进入被测器件(DUT)的散射电压和入射电压波之比，其反映了电路中各端口之间的信号关系。通常，S参数矩阵是根据射频多端口器件(DUT)的入射、检测端口次序排列的S参数得到的，如S_ij表示的是入射端口为j，检测端口为i。i＝j表示反射，i≠j表示传输。因此，对于一个n端口的网络，就有n×n个参数值。

宏模型是被测电子设备(DUT)的数学建模表示。S参数和Y、Z参数一样可以作为建模的频域输入数据。通过将测量或模拟得到的频域或时域数据投影到适当的基函数空间得到，这种投影或近似称为宏模型建模。常用有理函数作为基函数来进行被测电子设备(DUT)的分布式系统的宏模型建模，但它不能保证模型的因果性。特别是当延时以指数形式嵌入频域数据(S、Y或Z参数)中时，模型的阶数将达到无穷大，因此，具有有限阶的有理函数不能直接拟合含有延时的频域数据。

如果能够提前从频域数据(如S参数)中提取出延时项，再利用有理函数拟合不含延时的部分，将很大程度提高系统建模的精度和效率。目前普遍使用的传统延时提取方法是S.Grivet-Talocia于2006年提出的基于Gabor变换的时-频域分解方法，其流程见附图1、附图2以及附图3。这种传统方法是先将频域数据S(ξ)通过Gabor变换映射到时频域

其中窗函数为高斯窗：

W_ω,τ(ξ)＝W(ξ-ω)e^-jξτ (2.a)

其中τ在物理意义上是指时间或延时，积分的上下限是根据频域数据的带宽Ω进行计算。

然后再从时频数据获得具有瞬时功率特性的时域数据，即对

沿频域取平均值，将二维时频数据转为一维时域结果：

ε(τ)的局部峰值与延时项有对应关系，延时项可通过对ε(τ)提取局部峰值对应的时间点得到，所以只需要求出幅度的相对值即可。

从以上过程可以看出，最终提取的频域数据的延时项只需要用具有瞬时功率特性的一维时域数据，却需要先进行时频计算得到二维数据再计算得到一维数据，故计算中存在大量冗余信息。并且Gabor时频变换非常耗时，计算需要进行多次乘法和多次傅里叶逆变换，使得这种提取方法计算效率很低，从而严重影响了建模的效率。

发明内容

本发明的目的在于提供一种能够从频域数据(如S参数)中高精度地提取延时项的方法，有效提高射频电子设备(DUT)的宏模型建模的效率。

本发明第一方面提供了一种频域数据的延时项的提取方法，其特征在于，包括以下步骤：

(1)在不同频率下对射频n端口器件进行测量获取多组频域数据，每组频域数据对应一个n×n的参数矩阵

n为自然数；

(2)获取所述多组频域数据的各频率点{w₁ ⁰，w₂ ⁰，...，w_max ⁰}，根据各频率点{w₁ ⁰，w₂ ⁰，...，w_max ⁰}，从与各组频域数据分别对应的参数矩阵中提取具有延时项的元素k₁₁，构建集合U₁：{k₁₁(jw₁ ⁰)，k₁₁(jw₂ ⁰)，...，k₁₁(jw_max ⁰)}；

(3)对步骤(2)中的频率点{w₁ ⁰，w₂ ⁰，...，w_max ⁰}进行采样数据点数为自然数N的等间隔采样，得到新的频率点{w₁，w₂，...，w_N}，对所述集合U₁：{k₁₁(jw₁ ⁰)，k₁₁(jw₂ ⁰)，...，k₁₁(jw_max ⁰)}进行线性插值得到新的具有延时项的元素集合V：{k₁₁(jw₁)，k₁₁(jw₂)，...，k₁₁(jw_N)}；

(4))对步骤(3)中的具有延时项的元素集合V：{k₁₁(jw₁)，k₁₁(jw₂)，...，k₁₁(jw_N)}进行快速自相关计算

R(α)＝Xcorr(k₁₁(jw_t)，t∈[1，N]

(5)将步骤(4)得到的R(α)与指数项

相乘计算P(α)：

其中，σ是高斯分布标准差，α是各频率点之差w_i-w_t，其中i，t∈[1，N]且

(6)对P(α)进行快速傅里叶逆变换IFFT，然后对IFFT的结果取模获得频域数据的平均时频能量ε(τ)，

所述平均时频能量ε(τ)的局部峰值对应的时间点即为集合U₁频域数据的延时项，

对每组频域数据对应的参数矩阵中的其它元素k₁₂...k_nn同理执行上述步骤(2)-(6)的处理，直至提取出所述多组频域数据的全部延时项。

在一种可实现方式中，在步骤(4)中，对具有延时项的元素集合V：{k₁₁(jw₁)，k₁₁(jw₂)，...，k₁₁(jw_N)}进行快速傅里叶变换FFT得到序列F，再对序列F取共轭得到序列F^*，将两序列F和F^*对应项相乘再进行快速逆傅里叶变换IFFT得到R(α)，即R(α)＝IFFT(FFT(K(jw))·conj(FFT(K(jw))))，其中conj表示取共轭复数。

在一种可实现方式中，在步骤(5)中，将α的每一项取平方并除以常数项-4σ²得到新的序列M：

再对M序列的每一项取指数得到指数项序列：

在一种可实现方式中，本发明所述的方法还包括：对获得的集合U₁频域数据的每个延时项的能量进行量化，其中，沿所述平均时频能量ε(τ)的局部峰值两边的局部谷值区间进行积分，得到频域数据中每个延时项对应的能量贡献：

其中，i表示第i个延时项，τ_i、τ_i+1表示第i个峰值两边的谷值点。

设定一个能量阈值δ，按能量大小筛选出主要延时项：

根据本发明的上述方法，所述频域数据为射频/微波电路中的网络参数S、Y、Z参数或地震波数据中的一种。

本发明的第二方面提供了一种频域数据的延时项的提取方法，其特征在于，包括以下步骤：

(1)在不同频率下对射频n端口器件进行测量获取多组频域数据，每组频域数据对应一个n×n的参数矩阵

n为自然数；

(2)获取所述多组频域数据的各频率点{w₁ ⁰，w₂ ⁰，...，w_max ⁰}，根据各频率点{w₁ ⁰，w₂ ⁰，...，w_max ⁰}，从与各组频域数据分别对应的参数矩阵中提取具有延时项的元素k₁₁，构建集合U₁：{k₁₁(jw₁ ⁰)，k₁₁(jw₂ ⁰)，...，k₁₁(jw_max ⁰)}；

(3)对步骤(2)中的频率点{w₁ ⁰，w₂ ⁰，...，w_max ⁰}进行采样数据点数为自然数N的等间隔采样，得到新的频率点{w₁，w₂，...，w_N}，对所述集合U₁：{k₁₁(jw₁ ⁰)，k₁₁(jw₂ ⁰)，...，k₁₁(jw_max ⁰)}进行线性插值得到新的具有延时项的元素集合V：{k₁₁(jw₁)，k₁₁(jw₂)，...，k₁₁(jw_N)}；

(4)对步骤(3)中的具有延时项的元素集合V：{k₁₁(jw₁)，k₁₁(jw₂)，...，k₁₁(jw_N)}进行快速卷积计算：

(5)将步骤(4)得到的R(α)与指数项

相乘计算P(α)：

其中，σ是高斯分布标准差，α是各频率点之差w_i-w_t，其中i，t∈[1，N]且

(6)对P(α)进行快速傅里叶逆变换IFFT，然后对IFFT的结果取模获得频域数据的平均时频能量ε(τ)，

所述平均时频能量ε(τ)的局部峰值对应的时间点即为集合U₁频域数据的延时项，

对每组频域数据对应的参数矩阵中的其它元素k₁₂...k_nn同理执行上述步骤(2)-(6)的处理，直至提取出所述多组频域数据的全部延时项。

优选的，在本发明第二方面的方法中，进行快速卷积计算包括：对具有延时项的元素集合V：{k₁₁(jw₁)，k₁₁(jw₂)，...，k₁₁(jw_N)}取共轭对称得到

再对{k₁₁(jw₁)，k₁₁(jw₂)，...，k₁₁(jw_N)}和

分别进行快速傅里叶变换FFT，将得到的两个序列对应项相乘后再进行快速逆傅里叶变换IFFT得到R(α)，即R(α)＝IFFT(FFT(K(jw))·FFT(K^*(-jw)))。

优选的，在本发明第二方面的方法中，在步骤(5)中，

将α的每一项取平方并除以常数项-4σ²得到新的序列M：

再对M序列的每一项取指数得到指数项序列：

优选的，本发明的上述方法还包括：对获得的集合U1频域数据的每个延时项的能量进行量化，其中，沿所述平均时频能量ε(τ)的局部峰值两边的局部谷值区间进行积分，得到频域数据中每个延时项对应的能量贡献：

其中，i表示第i个延时项，τ_i、τ_i+1表示第i个峰值两边的谷值点。

设定一个能量阈值δ，按能量大小筛选出主要延时项：

本发明的第三方面还提供一种频域数据的延时项的提取装置，其特征在于，包括存储器和处理器，所述存储器上储存有在所述处理器上运行的程序，所述处理器运行所述程序时执行上述的频域数据的延时项的提取方法。

本发明第第四方面还提供一种计算机可读存储介质，其上存储有计算机指令，其特征在于，所述计算机指令被处理器运行时执行上述的频域数据的延时项的提取方法。

本发明利用近似方法巧妙避开传统计算方法中耗时较大、冗余信息较多的Gabor时频变换。根据S参数数据的带宽较宽，将原计算公式中的频率积分上下限近似为无穷，再运用高斯分布的概率密度积分公式对原计算式进行化简，直接使用输入的频域数据得到具有瞬时功率特性的时域结果。利用本发明所述方法提取得到的延时与传统Gabor变换方法得到的延时精确度相当，但在很大程度上提升了计算效率，大大节省了运算成本。

附图说明

图1是传统延时提取方法流程图。

图2是运用本发明方法的延时提取流程图。

图3是传统Gabor变换方法主要计算流程。

图4是本发明所述S参数延时项的提取方法主要计算流程图。

图5是快速卷积计算流程。

图6是快速自相关计算流程。

图7是传统Gabor变换方法与本发明所述方法得到的相对平均能量曲线对比原图，点线为传统方法，实线为本发明所述方法。

图8是图7所述相对平均能量曲线对比原图的局部放大图，点线为传统方法，实线为本发明所述方法。

具体实施方式

下面将结合本申请实施例中的附图，对本申请实施例中的技术方案进行清楚地描述，显然，所描述的实施例是本申请一部分实施例，而不是全部的实施例。基于本申请中的实施例，本领域普通技术人员获得的所有其他实施例，都属于本申请保护的范围。

本申请的说明书和权利要求书中的术语“第一”、“第二”等是用于区别类似的对象，而不用于描述特定的顺序或先后次序。应该理解这样使用的数据在适当情况下可以互换，以便本申请的实施例能够以除了在这里图示或描述的那些以外的顺序实施，且“第一”、“第二”等所区分的对象通常为一类，并不限定对象的个数，例如第一对象可以是一个，也可以是多个。此外，说明书以及权利要求中“和/或”表示所连接对象的至少其中之一，字符“/”，一般表示前后关联对象是一种“或”的关系。

本发明的主要思想是直接使用输入的频域数据得到具有瞬时功率特性的时域结果，具体地，本发明使用频域数据S(jξ)直接计算平均时频能量ε(τ)，具体过程如下。

首先将公式(3)：

展开得到：

再将公式(1)和(2)代入(4)得：

令

则有：

将公式(5)代入最里层括号得到：

其中

根据高斯分布概率密度积分公式得到：

将公式(8)的积分上下限近似为(-∞,+∞)，再由公式(9)得到A近似为常数，则公式(6)可近似得到：

因为最终只需要求ε(τ)的相对值，所以忽略公式中的常数系数，可得：

令ξ₁-ξ₂＝α，ξ₂＝ξ₁-α则有：

ε(τ)≈-∫P(α)·e^jατdα (12)

其中

R(α)＝∫S(ξ₁)S^*(ξ₁-α)dξ₁ (14)

根据公式(12)，ε(τ)可通过对P(α)进行逆傅里叶变换得到：

R(α)可以通过自相关或卷积计算获得：

根据公式(13)、(15)、(16)，从频域数据S(ξ)求出最终的平均时频能量ε(τ)，只需要进行一次卷积或自相关计算和一次快速傅里叶逆变换(IFFT)。如附图4、附图5所示，卷积和自相关可以利用FFT和IFFT快速计算，这样可以进一步降低计算代价，使得运算速度更快。

除了本文所说的Y、Z、S参数外，本发明所述方法还可应用于对其他频域数据的延时提取，比如对地震波的频域数据进行延时提取也可以应用本发明方法。。

下面结合图4，以射频/微波电路中的网络参数S说明本发明提取频域数据的延时项的具体实施步骤。可以理解的是，下面所述的方法同样可用于射频/微波电路中的网络参数Y、Z参数或其他频域数据的延时提取。

本发明提取频域数据的延时项的步骤包括：

S1)通过矢量分析仪等在不同频率下对例如射频多端口器件(DUT)进行测量，获得多组S参数数据。优选的，可以将所述获得多组S参数数据存储为Touchstone文件等。

在另一种可优选的方式中，也可通过对射频多端口器件进行数值模拟获得S参数数据，或者通过与远程服务器通信等获得不同频率下的被测射频设备的多组S参数数据。

S2)获取基于步骤S1得到的多组S参数数据的各频率点{w₁ ⁰，w₂ ⁰，...，w_max ⁰}，并对每个频率点所对应的S参数数据的S参数矩阵提取具有延时项的元素(比如S₁₁)，将提取出的具有延时项的元素组成集合U：{s₁₁(jw₁ ⁰)，s₁₁(jw₂ ⁰)，...，s₁₁(jw_max ⁰)}；

S3)设采样数据点数为N，对S参数数据的频率点{w₁ ⁰，w₂ ⁰，...，w_max ⁰}进行等间隔采样得到新的频率点{w₁，w₂，...，w_N}；根据该新的频率点{w₁，w₂，...，w_N}，对具有延时项的元素的所述集合U：{s₁₁(jw₁ ⁰)，s₁₁(jw₂ ⁰)，...，s₁₁(jw_max ⁰)}进行线性插值，得到新的具有延时项的元素的集合V：{s₁₁(jw₁)，s₁₁(jw₂)，...，s₁₁(jw_N)}。

S4)对通过步骤(S3)获得的新的具有延时项的集合V：{s₁₁(jw₁)，s₁₁(jw₂)，...，s₁₁(jw_N)}进行快速自相关或快速卷积计算：

式中的S表示矩阵或向量。该计算过程使用FFT和IFFT的快速卷积和快速互相关计算得到，具体方法如附图5、附图6所示。

采用快速自相关计算，即R(α)＝Xcorr(S(jw))的具体步骤为：

对具有延时项的元素集合V：{s₁₁(jw₁)，s₁₁(jw₂)，...，s₁₁(jw_N)}进行快速傅里叶变换FFT得到序列F，再对序列F取共轭得到序列F^*，将两序列F和F^*对应项相乘再进行快速逆傅里叶变换IFFT得到R(α)，即R(α)＝IFFT(FFT(S(jw))·conj(FFT(S(jw))))，其中conj表示取共轭复数。

采用快速卷积计算，即

的具体步骤为：

对具有延时项的元素集合V：{s₁₁(jw₁)，s₁₁(jw₂)，...，s₁₁(jw_N)}取共轭对称得到

再对{s₁₁(jw₁)，s₁₁(jw₂)，...，s₁₁(jw_N)}和

分别进行快速傅里叶变换FFT，将得到的两个序列对应项相乘后再进行快速逆傅里叶变换IFFT得到R(α)，即R(α)＝IFFT(FFT(S(jw))·FFT(S^*(-jw)))。

S5)将步骤S4)得到的R(α)与指数项

相乘得到P(α)。

其中σ是自定义的高斯分布标准差，一般定义为较小的常数，α是各频率点之差w_i-w_t，其中i，t∈[1，N]且

即α＝{w₁-w_N，w₁-w_N-1，...，0，w_N-w_N-1，w_N-w_N-2，...，w_N-w₁}，将α的每一项取平方并除以常数项-4σ²得到新的序列

再对M序列的每一项取指数得到指数项序列：

S6)计算S参数的平均时频能量ε(τ)。

即对P(α)进行快速傅里叶逆变换(IFFT)，因为得到的IFFT应该为实数，故对IFFT的结果取模得到平均时频能量ε(τ)。

S7)最后求出平均时频能量ε(τ)的局部峰值对应的时间点即为集合U的延时项。

经过多次实验证明，本发明所述的方法在保证频域数据的延时项提取的准确性的前提下，计算量以及计算复杂度都比传统的基于Gabor变换的时-频域分解方法小，提高了频域数据的延时项提取算法的整体性能。

a)本发明所述的方法窗长近似为无限长，与传统Gabor变换延时提取方法相比运算更方便，不需要调整窗长等信息，并且计算量和运算复杂度也比传统方法少很多。

当已知数据长度，本发明所述的方法可直接进行计算，而传统的基于Gabor变换的时-频域分解方法还需要设定窗长以及窗滑动的距离(hopsize)等参数，不同的窗长和窗滑动距离对其计算的速度以及最终的分辨率有很大影响，需要在保证分辨率的情况下用最快的速度计算，传统基于Gabor变换的时-频域分解方法的参数如下：

表1传统方法计算参数

窗长	overlap	hopsize	帧数	数据长度
					W	W-1	P	M＝ceil(N/P)	N

其中ceil是指向上取整。本例使用快速卷积的方法计算公式(16)。

对长度为n的数据进行IFFT和长度为n的两个数据补零后进行附图5的快速卷积，上述步骤所包含的乘法和加法次数如下表所示：

表2 IFFT和卷积所含乘法、加法次数

当输入频域数据例如S参数数据长度为N时，结合参数表1和表2，两种方法总共所用乘法和加法次数对比如下：

表3传统基于Gabor变换的时-频域分解方法和本发明所述的方法计算量对比

为了和本发明方法进行更好的对比，两个方法使用相同的标准差σ，传统基于Gabor变换的时-频域分解方法的窗长设为W＝N，窗滑动距离设为P＝1，则其帧数为M＝N，再将两个方法进行运算量的对比得到如下表：

表4相似条件下传统方法和本发明所述的方法计算量对比

粗略计算两种方法的计算复杂度：

表5相似条件下传统的基于Gabor变换的时-频域分解方法和本发明所述的方法计算复杂度对比

从表5中可看出，传统的基于Gabor变换的时-频域分解方法的乘法和加法复杂度是本发明所述方法的N倍，当数据量N非常大时，传统的基于Gabor变换的时-频域分解方法的计算量将远远超出本发明所述的方法，运算速率也远远低于本发明所述的方法。

使用不同长度的S参数数据在Matlab上运行两个算法，对比两个算法的运算时间，结果如下表所示：

表6传统的基于Gabor变换的时-频域分解方法与本发明所述的方法计算用时对比

从表中可知，本发明方法的运算时间受数据长度的影响比传统方法小很多，且用时比较少。随着数据长度增加，传统方法的运算时间增加的速度也越来越快。所以，当数据长度比较长时，本发明方法的运算速率比传统方法快很多。

b)本发明所述的方法求得的频域数据的延时项与传统的基于Gabor变换的时-频域分解方法准确度相当。

使用本发明所述的方法和传统的基于Gabor变换的时-频域分解方法对同一个频域数据例如S参数数据求延时项，数据长度为4001，窗长与数据长度一致，即等效为窗长无穷。将两种方法得到的相对平均能量绘制在一起，该曲线具有瞬时功率特性，如图7，8所示：

局部峰值的位置就是相应的延时项，对得到的具有瞬时功率特性的结果沿每个峰值两边的谷值构成的时间区间进行积分，求出每个峰值的量化能量：

其中，Ω_k表示第k个峰值两边谷值之间的区域。

根据每个峰值区域能量对整体能量的贡献大小来判断主要延时项，将能量占比大于阈值δ的峰值对应的延时作为主要延时项：

阈值δ根据数据设定，一般选择较小的值，如0.01～0.05之间，被忽略的延时项能量贡献很小，不会对最终的模型精度造成明显的影响。

此外，得到的主要延时项还可以进行后续的优化处理进而得出更加精确的延时，比如使用一些非线性优化方法。

运用上述方法，传统的基于Gabor变换的时-频域分解方法和本发明方法得到的局部峰值对应能量占比如下表所示：

表7四个主要峰值区域能量占比对比

本例选择阈值δ为0.01，则根据表7可得到3个主要的延时项，其时间位置如下表所示：

表8传统方法与本发明所述的方法求得主要峰值的位置

根据表8可以看出本发明所述的方法和传统的基于Gabor变换的时-频域分解方法求得的主要延时项误差不超过1.2％，这个误差在可接受范围内。

本发明还提供一种频域数据的延时项的提取装置，所述装置可以通过一般的计算机或计算机系统来实现，所述计算机或计算机系统具有输入装置、显示装置、外部I/F、通信I/F、处理器以及存储器。这些各硬件以能够经由总线进行通信的方式彼此连接。装置，所述装置可以通过一般的计算机或计算机系统来实现，所述计算机或计算机系统具有输入装置、显示装置、外部I/F、通信I/F、处理器以及存储器。这些各硬件以能够经由总线进行通信的方式彼此连接。

其中，所述处理器为计算机内部CPU，或专用CPU、DSP等处理部件。所述存储器装置，包括计算机可读存储介质，如计算机只读存储器(Read-Only Memory，ROM)、随机存取存储器(Random Access Memory，RAM)、磁碟或者光盘等。

所述存储器上储存有在所述处理器上运行的程序，所述处理器运行所述程序时执行上述频域数据的延时项的提取方法。

本发明还提供一种可读存储介质，所述可读存储介质上存储有程序或指令，该程序或指令被处理器执行时实现上述频域数据的延时项的提取方法的各个过程，且能达到相同的技术效果，为避免重复，这里不再赘述。

上面结合附图对本申请的实施例进行了描述，但是本申请并不局限于上述的具体实施方式，上述的具体实施方式仅仅是示意性的，而不是限制性的，本领域的普通技术人员在本申请的启示下，在不脱离本申请宗旨和权利要求所保护的范围情况下，还可做出很多形式，均属于本申请的保护之内。

Claims (12)

1.一种频域数据的延时项的提取方法，其特征在于，包括以下步骤：

(1)在不同频率下对射频n端口器件进行测量获取多组频域数据，每组频域数据对应一个n×n的参数矩阵

n为自然数；

(2)获取所述多组频域数据的各频率点{w₁ ⁰，w₂ ⁰，...，w_max ⁰}，根据各频率点{w₁ ⁰，w₂ ⁰，...，w_max ⁰}，从与各组频域数据分别对应的参数矩阵中提取具有延时项的元素k₁₁，构建集合U₁：{k₁₁(jw₁ ⁰)，k₁₁(jw₂ ⁰)，...，k₁₁(jw_max ⁰)}；

(3)对步骤(2)中的频率点{w₁ ⁰，w₂ ⁰，...，w_max ⁰}进行采样数据点数为自然数N的等间隔采样，得到新的频率点{w₁，w₂，...，w_N}，对所述集合U₁：{k₁₁(jw₁ ⁰)，k₁₁(jw₂ ⁰)，...，k₁₁(jw_max ⁰)}进行线性插值得到新的具有延时项的元素集合V：{k₁₁(jw₁)，k₁₁(jw₂)，...，k₁₁(jw_N)}；

(4)对步骤(3)中的具有延时项的元素集合V：{k₁₁(jw₁)，k₁₁(jw₂)，...，k₁₁(jw_N)}进行快速自相关计算：

R(α)＝Xcorr(k₁₁(jw_t)，t∈[1，N]

(5)将步骤(4)得到的R(α)与指数项

相乘计算P(α)：

其中，σ是高斯分布标准差，α是各频率点之差w_i-w_t，其中i，t∈[1，N]且

(6)对P(α)进行快速傅里叶逆变换IFFT，然后对IFFT的结果取模获得集合U₁频域数据的平均时频能量ε(τ)，

所述平均时频能量ε(τ)的局部峰值对应的时间点即为集合U₁频域数据的延时项，

对每组频域数据对应的参数矩阵中的其它元素k₁₂...knn同理执行上述步骤(2)-(6)的处理，直至提取出所述多组频域数据的全部延时项。

2.根据权利要求1所述的频域数据的延时项的提取方法，其特征在于，在步骤(4)中，对具有延时项的元素集合V：{k₁₁(jw₁)，k₁₁(jw₂)，...，k₁₁(jw_N)}进行快速傅里叶变换FFT得到序列F，再对序列F取共轭得到序列F^*，将两序列F和F^*对应项相乘再进行快速逆傅里叶变换IFFT得到R(α)，即R(α)＝IFFT(FFT(K(jw))·conj(FFT(K(jw))))，其中conj表示取共轭复数。

3.根据权利要求1所述的频域数据的延时项的提取方法，其特征在于，在步骤(5)中，

将α的每一项取平方并除以常数项-4σ²得到新的序列M：

再对M序列的每一项取指数得到指数项序列：

4.根据权利要求1-3中任一项所述的频域数据的延时项的提取方法，其特征在于，还包括：对获得的集合U₁频域数据的每个延时项的能量进行量化，其中，沿所述平均时频能量ε(τ)的局部峰值两边的局部谷值区间进行积分，得到集合U₁频域数据中每个延时项对应的能量贡献：

其中，i表示第i个延时项，τ_i、τ_i+1表示第i个峰值两边的谷值点。

设定一个能量阈值δ，按能量大小筛选出主要延时项：

5.根据权利要求1所述的频域的延时项的提取方法，其特征在于，所述频域数据为射频/微波电路中的网络参数S、Y、Z参数或地震波数据中的一种。

6.一种频域数据的延时项的提取方法，其特征在于，包括以下步骤：

(1)在不同频率下对射频n端口器件进行测量获取多组频域数据，每组频域数据对应一个n×n的参数矩阵

n为自然数；

(2)获取所述多组频域数据的各频率点{w₁ ⁰，w₂ ⁰，...，w_max ⁰}，根据各频率点{w₁ ⁰，w₂ ⁰，...，w_max ⁰}，从与各组频域数据分别对应的参数矩阵中提取具有延时项的元素k₁₁，构建集合U₁：{k₁₁(jw₁ ⁰)，k₁₁(jw₂ ⁰)，...，k₁₁(jw_max ⁰)}；

(3)对步骤(2)中的频率点(w₁ ⁰，w₂ ⁰，...，w_max ⁰}进行采样数据点数为自然数N的等间隔采样，得到新的频率点{w₁，w₂，...，w_N}，对所述集合U₁：{k₁₁(jw₁ ⁰)，k₁₁(jw₂ ⁰)，...，k₁₁(jw_max ⁰)}进行线性插值得到新的具有延时项的元素集合V：{k₁₁(jw₁)，k₁₁(jw₂)，...，k₁₁(jw_N)}；

(4)对步骤(3)中的具有延时项的元素集合V：{k₁₁(jw₁)，k₁₁(jw₂)，...，k₁₁(jw_N)}进行快速卷积计算：

(5)将步骤(4)得到的R(α)与指数项

相乘计算P(α)：

其中，σ是高斯分布标准差，α是各频率点之差w_i-w_t，其中i，t∈[1，N]且

(6)对P(α)进行快速傅里叶逆变换IFFT，然后对IFFT的结果取模获得频域数据的平均时频能量ε(τ)，

所述平均时频能量ε(τ)的局部峰值对应的时间点即为集合U₁频域数据的延时项，

对每组频域数据对应的参数矩阵中的其它元素k₁₂...knn同理执行上述步骤(2)-(6)的处理，直至提取出所述多组频域数据的全部延时项。

7.根据权利要求6所述的频域数据的延时项的提取方法，其特征在于，对步骤(3)中的新的具有延时项的元素集合V：{k₁₁(jw₁)，k₁₁(jw₂)，...，k₁₁(jw_N)}进行快速卷积计算包括：对具有延时项的元素集合V：{k₁₁(jw₁)，k₁₁(jw₂)，...，k₁₁(jw_N)}取共轭对称得到

再对{k₁₁(jw₁)，k₁₁(jw₂)，...，k₁₁(jw_N)}和

分别进行快速傅里叶变换FFT，将得到的两个序列对应项相乘后再进行快速逆傅里叶变换IFFT得到R(α)，即R(α)＝IFFT(FFT(K(jw))·FFT(K^*(-jw)))。

8.根据权利要求6所述的频域数据的延时项的提取方法，其特征在于，在步骤(5)中，

将α的每一项取平方并除以常数项-4σ²得到新的序列M：

再对M序列的每一项取指数得到指数项序列：

9.根据权利要求6至8中任一项所述的频域数据的延时项的提取方法，其特征在于，还包括：对获得的集合U1频域数据的每个延时项的能量进行量化，其中，沿所述平均时频能量ε(τ)的局部峰值两边的局部谷值区间进行积分，得到频域数据中每个延时项对应的能量贡献：

其中，i表示第i个延时项，τ_i、τ_i+1表示第i个峰值两边的谷值点。

设定一个能量阈值δ，按能量大小筛选出主要延时项：

10.根据权利要求6所述的频域的延时项的提取方法，其特征在于，所述频域数据为射频/微波电路中的网络参数S、Y、Z参数或地震波数据中的一种。

11.一种频域数据的延时项的提取装置，其特征在于，包括存储器和处理器，所述存储器上储存有在所述处理器上运行的程序，所述处理器运行所述程序时执行权利要求1-10任一项所述的频域数据的延时项的提取方法。

12.一种计算机可读存储介质，其上存储有计算机指令，其特征在于，所述计算机指令被处理器运行时执行权利要求1-10任一项所述的频域数据的延时项的提取方法。

Images