CN116953784A - 一种冰层振源信号到时拾取方法和装置 - Google Patents

Abstract

本发明涉及一种冰层振源信号到时拾取方法和装置，尤其涉及一种基于熵函数分析的到时拾取方法，属于极地声学技术领域。现有的到时提取方法很难满足冰层振源信号的精准到时拾取需求。本发明对采集的冰声信号进行分段处理分析，对每段信号依次进行模糊熵分析，截取模糊熵值突增的信号段；针对极地冰层信号，设计频域敏感的AIC算法特征函数，并对截取信号进行分析，以准确拾取振源信号到时。在没有丢失有效信号的基础上，更为精准的估计了冰层有效信号的到达临界时窗，能有效提高后续到时拾取以及信号特征分析的准确性，有效提高了冰层振源信号到时拾取的精度。

Description

一种冰层振源信号到时拾取方法和装置

技术领域

本发明涉及一种冰层振源信号到时拾取方法和装置，尤其涉及一种基于熵函数分析的到时拾取方法，属于极地声学技术领域。

背景技术

信号到时拾取是信号分析处理的基础，其精确程度直接影响了振源定位、介质波速测量的准确性。受冰层波导影响，冰层接收信号呈现明显的非平稳以及频散特性，信号峰值对应的频率、波速波动较大，信号相关法同样很难获取准确到时差。对于非平稳信号，到时拾取通常由具有丰富经验的人进行拾取，但其拾取速度和精度受人为因素影响较大，随着数据量增大，人工拾取已无法满足大量振动信号拾取的实时性要求。

目前国内外针对非平稳振动信号到时拾取的方法主要有长短时窗平均能量比法(STA/LTA)、赤池信息量准则(AIC)算法等。但对于冰层接收信号，其初始到达信号幅值较小，与噪声之间不存在幅值突变，因此基于能量变化的STA/LTA方法的拾取精度较差；而AIC方法需提前确定信号初始到达的时窗位置和时窗长度，否则到时拾取误差很大。专利申请《一种基于长短时窗和AR模型方差激增效应的P波到时拾取方法》提出使用STA/LTA方法初步确定信号到达时窗位置，通过截取信号片段结合AIC算法实现精确P波到时拾取。但受极地环境下的复杂噪声干扰影响，初始到达信号信噪比较低，导致STA/LTA方法很难通过选取阈值去有效确定信号初步到达时窗，同样降低了到时拾取的准确性。针对信噪比较低的信号，专利申请《一种微震震动波到时拾取方法》以及《一种岩体试件超声波初至拾取方法》中提出的先滤波后拾取到时信息的方法，但由于冰层接收信号与极区噪声信号之间还存在时频域混叠，会导致在滤波时丢失有效初至信号，增加到时拾取误差，因此也不适用。

综上所述，现有的到时提取方法很难满足冰层振源信号的精准到时拾取需求。因此，需针对极地冰声信号特点提出一种冰层振源信号到时拾取方法，以解决极地声源定位以及冰层波速测量的需求。

发明内容

本发明的目的是为了突破现有振源信号到时拾取方法面对冰层振源信号时适用性受限，无法有效提取有效信号的到时问题而提出的一种冰层振源信号到时拾取方法。该方法基于模糊熵对极地冰声信号复杂度进行分析，从而有效区分冰层噪声信号和有效信号，实现极地有效信号的到时拾取。

本发明对采集的冰声信号进行分段处理分析，对每段信号依次进行模糊熵分析，截取模糊熵值突增的信号段；针对极地冰层信号，设计频域敏感的AIC算法特征函数，并对截取信号进行分析，以准确拾取振源信号到时。

本发明提出一种冰层振源信号到时拾取方法，该方法包括以下步骤：

S1、对冰面上的接收信号X进行分段，得到m段长度为n的分段信号s_i,i＝1,2,…,m；

S2、计算每个分段信号s_i的模糊熵值FuEn_i；

S3、初至时窗截取

对每一分段信号s_i的模糊熵值FuEn_i做差求得模糊熵变化率v_i，通过最大值搜寻得到模糊熵变化率v_i的最大值v_j，将最大值v_j对应的分段信号以及前两段分段信号s_j-2,s_j-1,s_j组合得到信号长度为3n的初至振源信号段sig＝[s_j-2,s_j-1,s_j]；

S4、AIC函数值计算

设计AIC特征函数CF：

CF(sig(p))＝2sig(p)-(sig(p+1)+sig(p-1))，2＜p＜3n-1

计算AIC函数值Aic：

Aic(q)＝q*log{var(CF(sig([1,q])))}

+(3n-q-1)*log{var(CF(sig([q+1,3n])))}

其中q＝2,3,4,…,3n-2，var表示方差；

S5、振源信号到时拾取

记AIC函数值Aic最小值对应信号时刻为t，则t为振源信号的初至到达拾取时刻。

有利地，S1中，对信号长度为N的接收信号做傅里叶变换，确定有效信号的中心频率f₀，结合采样率fs，以10倍中心频率f₀的周期长度对信号X进行分段，得到分段信号s_i。

有利地，长度n＝10fs/f₀。

有利地，S2包括以下步骤：

S2.1、相空间重组

将每个分段信号s_i以相空间维数M进行相空间重组得到k个时间序列向量其中k＝1,2,3,…,n-M+1；

S2.2、计算最大绝对距离以及重组向量模糊隶属值

根据下式计算所有两组时间序列向量和/>之间的最大绝对距离/>k,l＝1,2,3，…,n-M+1,k≠l：

根据下式计算所有两组时间序列向量和/>之间的模糊隶属值/>

r为相似度容限；

S2.3、模糊度计算

针对每个时间序列向量求其模糊隶属值/>的特定平均值/>

通过下式计算分段信号s_i在相空间维数M下的模糊度

S2.4、计算模糊熵

令相空间维数M＝M+1，重复S2.1～S2.3得到分段信号s_i在相空间维数M+1下的模糊度计算每一分段信号s_i的模糊熵值：

有利地，相空间维数M为根据经验的给定值。

有利地，r＝0.1*σ，σ为分段信号s_i的标准差。

有利地，S3中，

本发明还提供一种冰层振源信号到时拾取装置，该装置包括检波器以及信号处理装置，通过检波器在冰面上获取接收信号X，通过信号处理装置对接收信号X进行处理和分析，执行如权利要求1-7中任一项所述的冰层振源信号到时拾取方法。

有利地，所述信号处理装置包括分段模块、模糊熵值计算模块、初至时窗截取模块、AIC值计算模块以及振源信号到时拾取模块，其中

所述分段模块对冰面上的接收信号X进行分段，得到m段长度为n的分段信号s_i,i＝1,2,…,m；

所述模糊熵值计算模块计算每个分段信号s_i的模糊熵值FuEn_i；

所述初至时窗截取模块对每一分段信号s_i的模糊熵值FuEn_i做差求得模糊熵变化率v_i，通过最大值搜寻得到模糊熵变化率v_i的最大值v_j，将最大值v_j对应的分段信号以及前两段分段信号s_j-2,s_j-1,s_j组合得到信号长度为3n的初至振源信号段sig＝[s_j-2,s_j-1,s_j]；

所述AIC值计算模块计算AIC特征函数CF：

CF(sig(p))＝2sig(p)-(sig(p+1)+sig(p-1))，2＜p＜3n-1

计算AIC函数值Aic：

Aic(q)＝q*log{var(CF(sig([1,q])))}

+(3n-q-1)*log{var(CF(sig([q+1,3n])))}

其中q＝2,3,4,…,3n-2，var表示方差；

所述振源信号到时拾取模块记AIC函数值Aic最小值对应信号时刻为t，则t为振源信号的初至到达时刻。

有利地，所述模糊熵值计算模块执行以下步骤：

S1、相空间重组

将每个分段信号s_i以相空间维数M进行相空间重组得到k个时间序列向量其中k＝1,2,3,…,n-M+1；

S2、计算最大绝对距离以及重组向量模糊隶属值

根据下式计算所有两组时间序列向量和/>之间的最大绝对距离/>k,l＝1,2,3，…,n-M+1,k≠l：

根据下式计算所有两组时间序列向量和/>之间的模糊隶属值/>

r为相似度容限；

S3、模糊度计算

针对每个时间序列向量求其模糊隶属值/>的特定平均值/>

通过下式计算分段信号s_i在相空间维数M下的模糊度

S4、计算模糊熵

令相空间维数M＝M+1，重复S2.1～S2.3得到分段信号s_i在相空间维数M+1下的模糊度计算每一分段信号s_i的模糊熵值：

与现有技术相比，本发明的有益效果是：

1.本发明基于冰层振源信号频域覆盖范围广的特点，利用模糊熵对信号混乱程度进行评估，在没有丢失有效信号的基础上，更为精准的估计了冰层有效信号的到达临界时窗，能有效提高后续到时拾取以及信号特征分析的准确性。

2.本发明针对初至冰层振源信号信噪比低、频域广特点，设计了对频率变化更为敏感、对幅值变化感知较弱的特征函数计算赤池信息量准则，有效提高了冰层振源信号到时拾取的精度。

已经讨论的特征、功能和优点可在各种示例中独立实现，或者可以在其他示例中进行组合。

附图说明

当结合附图阅读时，通过参考以下对本发明示例的详细描述，将最好地理解示例以及优选的使用模式、其他目的机器描述，其中：

图1是本发明方法是实施流程框图；

图2是本发明实例测试中实测振源信号时域波形图；

图3是本发明实例测试中实测振源信号局部放大时域波形图；

图4是本发明实例测试中信号分段及对应的模糊熵值图；

图5是本发明实例测试中模糊熵值的变化率图；

图6是本发明实例测试中振源到时拾取结果图；

图7是本发明实例对比中未使用本发明方法的振源到时拾取结果图。

具体实施方式

将参照附图更充分地描述所公开的示例，在附图中示出了所公开示例中的一些(但并非全部)。事实上，可描述许多不同的示例并且这些示例不应被解释为限于本文中阐述的示例。相反，描述这些示例，是的本发明公开彻底且完全，并且将把本发明公开的范围充分传达给本领域技术人员。

下面结合附图与具体实施方式对本发明作进一步详细描述。

图1所示为本发明方法的处理流程图。

步骤1、对冰面上的接收信号X进行分段

将检波器放置在冰面上获取接收信号X＝[x₁,x₂,x₃,…,x_N]，信号长度为N，对接收信号做傅里叶变换，可确定有效信号的中心频率f₀，结合采样率fs，以10倍中心频率f₀的周期长度对信号X进行分段，得到分段信号s_i。

操作时，将接收信号X以长度n(n＝10fs/f₀)均分为m段，得到m段长度为n的分段信号s_i,i＝1,2,…,m。

模糊熵是对信号复杂度的描述，在振源信号出现时，会出现明显熵增的特点，因此可以根据模糊熵值随分段信号的变化率，截取变化最大的信号段及前序分段信号，作为信号初至时窗。

步骤2、分段信号s_i的模糊熵值计算

2.1、相空间重组

将每个分段信号s_i以相空间维数M(给定值)进行相空间重组得到k个时间序列向量，其中k＝1,2,3,…,n-M+1。

2.2、计算最大绝对距离以及重组向量模糊隶属值

定义任意两组时间序列向量和/>(k,l＝1,2,3,…,n-M+1,k≠l)之间的最大绝对距离为/>并引入模糊隶属函数

其中r为相似度容限，通常取值为r＝0.1*σ，σ为分段信号s_i的标准差，计算所有两组时间序列向量组合的模糊隶属函数值。

2.3、计算模糊度

针对每个时间序列向量，求其模糊隶属值的特定平均值/>

定义为分段信号s_i在相空间维数M下的模糊度

计算每个分段信号s_i在相空间维数M下的模糊度

2.4、计算模糊熵

令相空间维数M＝M+1，重复步骤2.1～2.3得到分段信号s_i在相空间维数M+1下的模糊度计算每一分段信号s_i的模糊熵值：

步骤3、初至时窗截取

对每一分段信号s_i的模糊熵值FuEn_i做差求得模糊熵变化率v_i，即

通过最大值搜寻得到模糊熵变化率v_i的最大值v_j，根据模糊熵值突变结果可知分段信号中包含了振源信号，将最大值v_j对应的分段信号以及前两段分段信号s_j-2,s_j-1,s_j组合得到信号长度为3n的初至振源信号段sig＝[s_j-2,s_j-1,s_j]；

为更为精准的提取信号到时，本发明应用赤池信息量准则(AIC)算法进行计算。针对初至信号信噪比低的问题，本发明针对冰层振源信号信噪比低，频率变化更为明显的特点，设计了频率变化敏感度高，幅值敏感度低的特征函数，以此带入AIC函数，依据AIC函数最小值即信号到达时刻，求取准确的信号到时。

步骤4、AiC函数值计算

设计AIC特征函数CF：

CF(sig(p))＝2sig(p)-(sig(p+1)+sig(p-1))，2＜p＜3n-1

计算AIC函数值Aic：

Aic(q)＝q*log{var(CF(sig([1,q])))}

+(3n-q-1)*log{var(CF(sig([q+1,3n])))}

其中q＝2,3,4,…,3n-2，var表示方差。

步骤5、振源信号到时拾取

记AIC函数值Aic最小值对应信号时刻为t，则t为振源信号的初至到达时刻。

下面以松花江实测信号进一步对本发明方法的可行性进行说明。

如图2所示为松花江实测的震源接收信号，可以看到噪声整体频域比信号要低。放大信号初至时段如图3所示，接收信号频率刚变化时信噪比很低，人工选取到时为365.6ms。对接收信号进行分段并求取模糊熵，其结果如图4所示，可以看到分段后，在噪声区域的模糊熵值普遍较低，有效信号区域，模糊熵值都较高。求得模糊熵值的变化率如图5所示，模糊熵变化率峰值在第8段信号，根据上述步骤3选取第5～8段信号作为振源信号初至时窗，如图6中实线所示。通过本发明设计的特征函数，计算的AIC函数值如图6中虚线所示，AIC函数的最小值点对应的到时即为信号到达时间，本发明方法的到时估计时间为365.1ms，误差为0.5ms。以AIC函数默认的频域敏感特征函数计算结果如图7所示，到时估计为376.8ms，误差为11.2ms。因此本方法可以更为精确的估计冰层振源信号到达时间。

已出示示例和描述的目的展示了对不同有利布置的描述，但是该描述并不旨在是排他性的或限于所公开形式的示例。许多修改形式和变化形式对于本领域技术人员而言将是显而易见的。另外，不同的有利示例可描述与其他有利示例相比具有不同的有点。选择和描述所选择的一个示例或多个示例，以便最佳地说明示例的原理、实际应用，并且使本领域的普通技术人员能够理解本公开有进行了适于所料想特定使用的各种修改的各种示例。

Claims (10)

1.一种冰层振源信号到时拾取方法，其特征在于，该方法包括以下步骤：

S1、对冰面上的接收信号X进行分段，得到m段长度为n的分段信号s_i,i＝1,2,…,m；

S2、计算每个分段信号s_i的模糊熵值FuEn_i；

S3、初至时窗截取

对每一分段信号s_i的模糊熵值FuEn_i做差求得模糊熵变化率v_i，通过最大值搜寻得到模糊熵变化率v_i的最大值v_j，将最大值v_j对应的分段信号以及前两段分段信号s_j-2,s_j-1,s_j组合得到信号长度为3n的初至振源信号段sig＝[s_j-2,s_j-1,s_j]；

S4、AiC函数值计算

设计AIC特征函数CF：

CF(sig(p))＝2sig(p)-(sig(p+1)+sig(p-1))，2＜p＜3n-1

计算AIC函数值Aic：

Aic(q)＝q*log{var(CF(sig([1,q])))}

+(3n-q-1)*log{var(CF(sig([q+1,3n])))}

其中q＝2,3,4,…,3n-2，var表示方差；

S5、振源信号到时拾取

记AIC函数值Aic最小值对应信号时刻为t，则t为振源信号的初至到达拾取时刻。

2.根据权利要求1所述的冰层振源信号到时拾取方法，其特征在于：S1中，对信号长度为N的接收信号做傅里叶变换，确定有效信号的中心频率f₀，结合采样率fs，以10倍中心频率f₀的周期长度对信号X进行分段，得到分段信号s_i。

3.根据权利要求2所述的冰层振源信号到时拾取方法，其特征在于：长度n＝10fs/f₀。

4.根据权利要求1所述的冰层振源信号到时拾取方法，其特征在于：S2包括以下步骤：

S2.1、相空间重组

将每个分段信号s_i以相空间维数M进行相空间重组得到k个时间序列向量其中k＝1,2,3,…,n-M+1；

S2.2、计算最大绝对距离以及重组向量模糊隶属值

根据下式计算所有两组时间序列向量和/>之间的最大绝对距离/>

根据下式计算所有两组时间序列向量和/>之间的模糊隶属值/>

r为相似度容限；

S2.3、模糊度计算

针对每个时间序列向量求其模糊隶属值/>的特定平均值/>

通过下式计算分段信号s_i在相空间维数M下的模糊度

S2.4、计算模糊熵

令相空间维数M＝M+1，重复S2.1～S2.3得到分段信号s_i在相空间维数M+1下的模糊度计算每一分段信号s_i的模糊熵值：

5.根据权利要求4所述的冰层振源信号到时拾取方法，其特征在于：相空间维数M为根据经验的给定值。

6.根据权利要求4所述的冰层振源信号到时拾取方法，其特征在于：r＝0.1*σ，σ为分段信号s_i的标准差。

7.根据权利要求1所述的冰层振源信号到时拾取方法，其特征在于：S3中，

8.一种冰层振源信号到时拾取装置，其特征在于：该装置包括检波器以及信号处理装置，通过检波器在冰面上获取接收信号X，通过信号处理装置对接收信号X进行处理和分析，执行如权利要求1-7中任一项所述的冰层振源信号到时拾取方法。

9.根据权利要求8所述的冰层振源信号到时拾取装置，其特征在于：所述信号处理装置包括分段模块、模糊熵值计算模块、初至时窗截取模块、AiC函数值计算模块以及振源信号到时拾取模块，其中

所述分段模块对冰面上的接收信号X进行分段，得到m段长度为n的分段信号s_i,i＝1,2,…,m；

所述模糊熵值计算模块计算每个分段信号s_i的模糊熵值FuEn_i；

所述初至时窗截取模块对每一分段信号s_i的模糊熵值FuEn_i做差求得模糊熵变化率v_i，通过最大值搜寻得到模糊熵变化率v_i的最大值v_j，将最大值v_j对应的分段信号以及前两段分段信号s_j-2,s_j-1,s_j组合得到信号长度为3n的初至振源信号段sig＝[s_j-2,s_j-1,s_j]；

所述AIC值计算模块计算ARC特征函数CF：

CF(sig(p))＝2sig(p)-(sig(p+1)+sig(p-1))，2＜p＜3n-1

计算AIC函数值Aic：

Aic(q)＝q*log{var(CF(sig([1,q])))}

+(3n-q-1)*log{var(CF(sig([q+1,3n])))}

其中q＝2,3,4,…,3n-2，var表示方差；

所述振源信号到时拾取模块记AIC函数值Aic最小值对应信号时刻为t，则t为振源信号的初至到达时刻。

10.根据权利要求9所述的冰层振源信号到时拾取装置，其特征在于：所述模糊熵值计算模块执行以下步骤：

S1、相空间重组

将每个分段信号s_i以相空间维数M进行相空间重组得到k个时间序列向量其中k＝1,2,3,…,n-M+1；

S2、计算最大绝对距离以及重组向量模糊隶属值

根据下式计算所有两组时间序列向量和/>之间的最大绝对距离/>

根据下式计算所有两组时间序列向量和/>之间的模糊隶属值/>

r为相似度容限；

S3、模糊度计算

针对每个时间序列向量求其模糊隶属值/>的特定平均值/>

通过下式计算分段信号s_i在相空间维数M下的模糊度

S4、计算模糊熵

令相空间维数M＝M+1，重复S2.1～S2.3得到分段信号s_i在相空间维数M+1下的模糊度计算每一分段信号s_i的模糊熵值：