本發(fā)明涉及地球物理勘探的地震數(shù)據(jù)處理技術(shù)領(lǐng)域，尤其涉及一種含噪地震信號初至走時(shí)拾取方法及裝置。

背景技術(shù)：

地震波在有限區(qū)域的介質(zhì)中的傳播需要消耗時(shí)間，地震信號由震源(激發(fā)點(diǎn))傳播到接收(檢波點(diǎn))的初至?xí)r刻數(shù)據(jù)蘊(yùn)含著介質(zhì)內(nèi)部的物性特征(狀態(tài)、密度等)。地震記錄初至波(p波)走時(shí)拾取是地震層析正演計(jì)算和靜校正技術(shù)中一項(xiàng)基礎(chǔ)而又重要的工作，對震源定位，油藏描述，探明地質(zhì)缺陷分布等具有重要意義。初至波走時(shí)的手工拾取已不符合現(xiàn)今大規(guī)模地震勘探數(shù)據(jù)的需求。一般的計(jì)算機(jī)自動拾取方法能滿足于地震記錄背景較平靜且初至波形簡單的地震資料。但對于具有較強(qiáng)背景噪聲干擾，初至波形變化較大，各種波相互干擾的地震記錄，往往難以獲得準(zhǔn)確的初至?xí)r間。因此，在運(yùn)用自動拾取算法前對地震隨機(jī)噪聲進(jìn)行壓制處理，是提高地震初至波走時(shí)拾取精度的有效手段。

技術(shù)實(shí)現(xiàn)要素：

本發(fā)明的目的在于提供一種含噪地震信號初至走時(shí)拾取方法及裝置，以解決自然或人工地震資料中存在隨機(jī)噪聲干擾時(shí)難以獲得準(zhǔn)確的初至?xí)r間的問題，該方法可提高地震信號信噪比，提升地震初至波走時(shí)拾取精度。

為達(dá)到上述目的，一方面，一種含噪地震信號初至走時(shí)拾取方法，其步驟：

s1：采用自適應(yīng)噪聲完備集合經(jīng)驗(yàn)?zāi)B(tài)分解(complementaryensembleempiricalmodedecompositionwithadaptivenoise,ceemdan)方法對原始含噪地震信號s(t)分解處理得到一系列從高頻到低頻的有限個(gè)本征模態(tài)函數(shù)(intrinsicmodefunction,imf)分量及殘余分量；

s2：隨機(jī)噪聲往往分布在第一個(gè)或前幾個(gè)高頻imf分量中，采用有效指標(biāo)函數(shù)(effectiveindexfunction,eif)確定各階imf分量中有效信號與噪聲的能量分界點(diǎn)，判斷含噪較多的高頻imf分量；

s3：依據(jù)能量分界點(diǎn)對部分高頻imf分量進(jìn)行小波閾值重構(gòu)小波系數(shù)，并保持低頻imf分量及殘余分量不變；

s4：將小波閾值中去噪后的高頻imf分量和不做處理的低頻imf分量以及殘余分量累加重構(gòu)，即得到去噪后的地震信號x(t)：

式中，n為ceemdan分解得到的imf分量級數(shù)；m為需要進(jìn)行小波閾值去噪的高頻imf分量數(shù)目；imfi′(t)為進(jìn)行小波閾值去噪后的高頻imf分量；imfi(t)為不需要處理的低頻imf分量；rn(t)為ceemdan分解后的殘余分量；

s5：采用改進(jìn)的滑動時(shí)窗能量比法識別地震事件并初步確定初至波走時(shí)為t0；

s6：去噪后地震序列可截?cái)酁槎鄠€(gè)局部平穩(wěn)部分，其各部分均可表示為自回歸(autorecursive,ar)模型。當(dāng)ar模型階數(shù)固定時(shí)，采用aic信息準(zhǔn)則法(akaikeinformationcriterion,aic)計(jì)算ar模型的aic函數(shù)最小值，該最小值位置，可充當(dāng)兩段地震序列的分界點(diǎn)。在對應(yīng)t0時(shí)刻前后取一時(shí)窗，使用aic函數(shù)計(jì)算在該時(shí)窗范圍內(nèi)的局部極小值，即得到精確的初至波走時(shí)t；

s7：根據(jù)約束準(zhǔn)則對拾取的初至波走時(shí)進(jìn)行評價(jià)，對含異常初至?xí)r間信息的地震道進(jìn)行約束處理。

所述原始含噪地震信號的ceemdan分解步驟包括：

s1_1，定義算子ek(·)為產(chǎn)生第k個(gè)imf模態(tài)分量的操作過程；wⁱ(t)表示零均值單位方差的正負(fù)白噪聲；i＝1，2，...，i；系數(shù)βk可在每一個(gè)模態(tài)分解階段選擇合適的信噪比(signal-noiseratio，snr)；k＝1，2，...，k表示分解階段；

s1_2，向原始地震信號s(t)中添加零均值單位方差的正負(fù)白噪聲wⁱ(t)，得到處理后信號序列sⁱ(t)為：sⁱ(t)＝s(t)+β0wⁱ(t)；

s1_3，對sⁱ(t)進(jìn)行i次驗(yàn)?zāi)B(tài)分解(empiricalmodedecomposition，emd)處理，得到第一階模態(tài)分量imf1(t)為：

s1_4，對于k＝1，計(jì)算第一階模態(tài)分量后殘余分量r1(t)為：r1(t)＝s(t)-imf1(t)；

s1_5，對信號r1(t)+β1e1(w¹(t))進(jìn)行emd分解，得到第二階模態(tài)分量為：

s1_6，對其余每個(gè)階段，即k＝1，2，...，k，與步驟s1_4的計(jì)算過程一致，首先計(jì)算第k階殘余分量，再計(jì)算第k+1階模態(tài)分量，即有：

s1_7，重復(fù)以上步驟直至殘余分量成為一個(gè)單調(diào)函數(shù)或整個(gè)時(shí)間序列上極大值點(diǎn)數(shù)和極小值點(diǎn)數(shù)之和不大于兩個(gè)，則結(jié)束分解過程。

所述eif有效指標(biāo)函數(shù)定義如下：

式中，x(t)為去噪后的地震信號，n代表信號長度；imf(t)為分解出的各級模態(tài)分量；eif越小就表示imf(t)分量越近似于原始信號，eif的較大值被認(rèn)為是噪聲。

所述小波閾值去噪步驟包括：首先將能量分界點(diǎn)劃分出的含噪高頻imf分量信號轉(zhuǎn)換至小波域，求取小波系數(shù)wj：

wj＝w(imf(t))；

其中w(·)表示小波變換，j為小波分解層數(shù)。

設(shè)置一個(gè)閾值，小于該閾值的小波系數(shù)被認(rèn)為是由噪聲產(chǎn)生的，將其去除。閾值處理采用軟閾值函數(shù)

式中，sgn()為符號函數(shù)，wj為小波系數(shù)，λ為閾值，其中σ為噪聲標(biāo)準(zhǔn)差，n為信號長度。噪聲標(biāo)準(zhǔn)差median(|wj|)為小波多分辨級分解系數(shù)的中值。

最后對小波系數(shù)wj進(jìn)行逆變換，即得到小波閾值去噪后的高頻imf分量信號。

所述改進(jìn)的滑動時(shí)窗能量比法的原理是：地震記錄初至波來臨之前，大部分為噪聲信號，當(dāng)初至波來臨時(shí)初至?xí)r間前后時(shí)窗內(nèi)的地震能量特征有較大的差異，該時(shí)刻前后時(shí)窗能量比值a(t)達(dá)到最大值，因此，拾取a(t)最大值對應(yīng)的時(shí)間即可獲得地震信號的初至走時(shí)。在前、后時(shí)窗能量中加入穩(wěn)定因子aw，增強(qiáng)其穩(wěn)定性，改進(jìn)的滑動時(shí)窗能量比法的計(jì)算公式可表示為：

式中，x(t)為地震記錄振幅值，α為穩(wěn)定系數(shù)；為地震道相對能量；t0為前時(shí)窗起點(diǎn)，t1為前時(shí)窗終點(diǎn)或后時(shí)窗起點(diǎn)，t2為后時(shí)窗終點(diǎn)。

所述去噪后地震信號的aic值采用下式計(jì)算：

式中，k為有效地震信號和隨機(jī)噪聲分量的最佳分界點(diǎn)，σ1和σ2則分別表示隨機(jī)噪聲和有效地震信號分量的方差，m表示ar模型階數(shù)，n為信號長度，c為常數(shù)。

在地震波初至來臨時(shí)波形發(fā)生較大變化，隨機(jī)噪聲和地震信號具有較大的數(shù)學(xué)統(tǒng)計(jì)差異，兩者擬合度最差，可求得最小aic值。因此，計(jì)算地震波的aic極小值可用于地震初至走時(shí)的確定。

所述初至拾取異常的道采用如下步驟處理：

根據(jù)同一炮的相鄰兩個(gè)檢波點(diǎn)的初至?xí)r間不會有較大突變，相繼各道合理的拾取時(shí)間應(yīng)形成折射段或同相軸的先驗(yàn)知識，對初至走時(shí)進(jìn)行評價(jià)，然后構(gòu)造約束準(zhǔn)則完成異常道的初至拾取。計(jì)算相鄰道的初至?xí)r間差絕對值和的平均值式中m為單炮總道數(shù)，ti是每道的初至走時(shí)，對相鄰初至?xí)r間差分進(jìn)行組合分析，如果滿足條件(-1≤n≤9,1≤i≤m,1≤j≤m)，則判定該道初至?xí)r間的可信度為1，如果不滿足，則判定該道初至?xí)r間的可信度為0。確定了部分初至?xí)r間拾取可信度為1的道，對篩選后被認(rèn)為可信度為0的道，通過道間樣條插值或者線性插值進(jìn)行重新確定，即得到異常道糾正處理后的初至走時(shí)。其中0＜k≤2，針對信噪比高的地震資料，k參數(shù)值可取較大值，反之，k取較小值。

一種含噪地震信號初至走時(shí)拾取裝置，包括：噪聲壓制單元u1、初至拾取單元u2。

噪聲壓制單元，用于壓制原始地震信號中的隨機(jī)噪聲，是地震初至波走時(shí)準(zhǔn)確拾取的前提條件；

初至拾取單元，用于拾取地震波初至來臨的時(shí)間信息。

噪聲壓制單元與初至拾取單元通過電路連接。

所述噪聲壓制單元包括：ceemdan分解模塊m1、imf分量分析模塊m2、小波閾值去噪模塊m3、分量重構(gòu)模塊m4；

ceemdan分解模塊，用于對原始含噪地震信號進(jìn)行逐級分解得到一系列從高頻到低頻的有限個(gè)imf分量以及殘余分量；

imf分量分析模塊，采用有效的指標(biāo)函數(shù)eif來確定能量分界點(diǎn)劃分分解得到的有效信號或噪聲(大部分的高頻噪聲會出現(xiàn)在前幾個(gè)模態(tài)分量中)：

式中，n代表信號長度；imf(t)為原信號分解出的各級模態(tài)分量；eif越小就表示imf(t)分量越近似于原始信號，eif的較大值被認(rèn)為是噪聲。

小波閾值去噪模塊，將imf分量分析模塊中劃分出的含噪高頻imf分量再次進(jìn)行小波閾值去噪；

分量重構(gòu)模塊，用于將去噪后的高頻imf分量和低頻imf分量以及殘余量進(jìn)行累加重構(gòu)，即得到去噪后的地震信號。

ceemdan分解模塊、imf分量分析模塊、小波閾值去噪模塊和分量重構(gòu)模塊依次按順序通過電路連接。

所述初至拾取單元包括：改進(jìn)滑動時(shí)窗能量比法預(yù)拾取模塊m5、aic法再拾取模塊m6和初至評價(jià)與異常處理模塊m7；

改進(jìn)滑動時(shí)窗能量比法預(yù)拾取模塊，用于識別地震事件并大致確定波至?xí)r刻，地震記錄初至波來臨之前，大部分為噪聲信號，當(dāng)初至波來臨時(shí)初至?xí)r間前后時(shí)窗內(nèi)的地震能量特征有較大的差異，該時(shí)刻前后時(shí)窗能量比值a(t)達(dá)到最大值，因此，拾取a(t)最大值對應(yīng)的時(shí)間即可獲得地震信號的初至走時(shí)。在前、后時(shí)窗能量中加入穩(wěn)定因子aw，增強(qiáng)其穩(wěn)定性，改進(jìn)的滑動時(shí)窗能量比法的計(jì)算公式可表示為：

式中，x(t)為地震記錄振幅值，α為穩(wěn)定系數(shù)；為地震道相對能量；t0為前時(shí)窗起點(diǎn)，t1為前時(shí)窗終點(diǎn)或后時(shí)窗起點(diǎn)，t2為后時(shí)窗終點(diǎn)。

aic法再拾取模塊，用于初至?xí)r刻的精確拾取，在改進(jìn)滑動時(shí)窗能量比法預(yù)拾取模塊的基礎(chǔ)上，于預(yù)拾取走時(shí)對應(yīng)時(shí)刻前后取一時(shí)窗，使用aic信息準(zhǔn)則法計(jì)算地震序列ar模型在該時(shí)窗范圍內(nèi)的局部極小值，該極小值對應(yīng)時(shí)間即為精確的初至波走時(shí)，aic值計(jì)算公式如下式：

式中，k為有效地震信號和隨機(jī)噪聲分量的最佳分界點(diǎn)，σ1和σ2則分別表示隨機(jī)噪聲和有效地震信號分量的方差，m表示ar模型階數(shù)，n為信號長度，c為常數(shù)。

在地震波初至來臨時(shí)波形發(fā)生較大變化，隨機(jī)噪聲和地震信號具有較大的數(shù)學(xué)統(tǒng)計(jì)差異，兩者擬合度最差，可求得最小aic值。因此，計(jì)算地震波的aic極小值可用于地震初至走時(shí)的確定。

初至評價(jià)與異常處理模塊，用于處理無信號道或干擾道等異常情況出現(xiàn)時(shí)，計(jì)算機(jī)自動初至拾取方法出錯的情況。根據(jù)同一炮的相鄰兩個(gè)檢波點(diǎn)的初至?xí)r間不會有較大突變，相繼各道合理的拾取時(shí)間應(yīng)形成折射段或同相軸的先驗(yàn)知識，對初至進(jìn)行評價(jià)，然后構(gòu)造約束準(zhǔn)則完成異常道的初至拾取。計(jì)算相鄰道的初至?xí)r間差絕對值和的平均值式中m為單炮總道數(shù)，ti是每道的初至走時(shí)，對相鄰初至?xí)r間差分進(jìn)行組合分析，如果滿足條件(-1≤n≤9,1≤i≤m,1≤j≤m)，則判定該道初至?xí)r間的可信度為1，如果不滿足，則判定該道初至?xí)r間的可信度為0。確定了部分初至?xí)r間拾取可信度為1的道，對篩選后被認(rèn)為可信度為0的道，通過道間樣條插值或者線性插值進(jìn)行重新確定，即得到異常道糾正處理后的初至走時(shí)。其中0＜k≤2，針對信噪比高的地震資料，k參數(shù)值可取較大值，反之，k取較小值。

改進(jìn)滑動時(shí)窗能量比法預(yù)拾取模塊、aic法再拾取模塊和初至評價(jià)與異常處理模塊依次按順序通過電路連接。

本發(fā)明的優(yōu)點(diǎn)：依據(jù)先驗(yàn)知識對拾取的初至?xí)r刻進(jìn)行評價(jià)，將評價(jià)可信度高的道作為可信度低道的約束條件采用插值計(jì)算進(jìn)行二次拾取等一系列操作，免人工干擾，抗干擾能力強(qiáng)，自動剔除低信噪比道，保證初至拾取結(jié)果的可靠性和提高拾取速度。本發(fā)明可運(yùn)用于不同復(fù)雜度的地震資料，靜校正效果好，初至拾取質(zhì)量高，對實(shí)際生產(chǎn)、研究應(yīng)用具有較高意義。

附圖說明

為了更清楚地說明本申請實(shí)施例或現(xiàn)有技術(shù)中的技術(shù)方案，下面將對實(shí)施例或現(xiàn)有技術(shù)描述中所需要使用的附圖作簡單地介紹，顯而易見地，下面描述中的附圖僅僅是本申請中記載的一些實(shí)施例，對于本領(lǐng)域普通技術(shù)人員來講，在不付出創(chuàng)造性勞動性的前提下，還可以根據(jù)這些附圖獲得其他的附圖。

附圖1是一種含噪地震信號初至走時(shí)拾取方法實(shí)現(xiàn)流程圖；

附圖2是一種含噪地震信號初至走時(shí)拾取裝置的模塊結(jié)構(gòu)示意圖；

附圖3是一種含噪地震信號初至走時(shí)拾取方法實(shí)現(xiàn)示意圖；

附圖4是本申請一實(shí)施例的合成地震道集數(shù)據(jù)去噪效果對比圖；

附圖5是本申請一實(shí)施例的合成地震道集數(shù)據(jù)中各道去噪后信噪比變化對比圖；

附圖6是本申請一實(shí)施例的合成地震道集數(shù)據(jù)中各道去噪后初至拾取效果對比圖；

附圖7是本發(fā)明一實(shí)施例的四組實(shí)際單炮地震記錄；

附圖8是本申請一實(shí)施例的四組實(shí)際單炮地震記錄的初至拾取分布示意圖；

附圖9是本申請一實(shí)施例的四組實(shí)際單炮地震記錄的初至拾取精度roc曲線示意圖；

具體實(shí)施方式

下面將結(jié)合本發(fā)明實(shí)施例中的附圖，對本發(fā)明實(shí)施例中的技術(shù)方案進(jìn)行清楚、完整地描述，顯然，所描述的實(shí)施例僅僅是本發(fā)明一部分實(shí)施例，而不是全部的實(shí)施例。基于本發(fā)明中的實(shí)施例，本領(lǐng)域普通技術(shù)人員在沒有作出創(chuàng)造性勞動前提下所獲得的所有其他實(shí)施例，都屬于本發(fā)明保護(hù)的范圍。

本實(shí)施例是一種含噪地震信號初至走時(shí)拾取方法及裝置，見附圖1至附圖10所示。在本發(fā)明實(shí)施例中，首先，原始含噪地震信號通過ceemdan分解得到各級imf分量，采用eif有效指標(biāo)函數(shù)確定有效信號和噪聲的能量分界點(diǎn)，對部分高頻imf分量再次進(jìn)行小波閾值去噪，并保持低頻imf分量和殘余分量不變；將去噪后的高頻imf分量和低頻imf分量以及殘余量進(jìn)行累加重構(gòu)，即得到去噪后的地震信號。然后，采用改進(jìn)滑動時(shí)窗能量比法在地震信號上取一滑動時(shí)窗計(jì)算前后時(shí)窗能量比值，當(dāng)能量比值大于設(shè)定閾值時(shí)，即可識別地震事件并大致確定波至?xí)r刻；在該時(shí)刻前后取一時(shí)窗，使用aic方法計(jì)算ar模型在該時(shí)窗范圍內(nèi)的局部極小值，即得到精確的初至波走時(shí)。最后，循環(huán)道集判斷是否有初至?xí)r間拾取異常的道出現(xiàn)，如有，采用異常處理措施進(jìn)行糾正處理。

如附圖1，本請實(shí)施例提供一種含噪地震信號初至走時(shí)拾取方法。該方法包括：

s1：采用自適應(yīng)噪聲完備集合經(jīng)驗(yàn)?zāi)B(tài)分解ceemdan方法對原始含噪地震信號s(t)分解處理得到一系列從高頻到低頻有限個(gè)的本征模態(tài)函數(shù)imf分量及殘余分量。原始含噪地震信號的ceemdan分解步驟包括：

s1_1，定義算子ek(·)為產(chǎn)生第k個(gè)imf模態(tài)分量的操作過程；wⁱ(t)表示零均值單位方差的正負(fù)白噪聲；i＝1，2，...，i；系數(shù)βk可在每一個(gè)模態(tài)分解階段選擇合適的信噪比snr；k＝1，2，...，k表示分解階段；

s1_2，向原始地震信號s(t)中添加零均值單位方差的正負(fù)白噪聲wⁱ(t)，得到處理后信號序列sⁱ(t)為：sⁱ(t)＝s(t)+β0wⁱ(t)；

s1_3，對sⁱ(t)進(jìn)行i次驗(yàn)?zāi)B(tài)分解emd處理，得到第一階模態(tài)分量imf1(t)為：

s1_4，對于k＝1，計(jì)算第一階模態(tài)分量后殘余分量r1(t)為：r1(t)＝s(t)-imf1(t)；

s1_5，對信號r1(t)+β1e1(w¹(t))進(jìn)行emd分解，得到第二階模態(tài)分量為：

s1_6，對其余每個(gè)階段，即k＝1，2，...，k，與步驟s1_4的計(jì)算過程一致，首先計(jì)算第k階殘余分量，再計(jì)算第k+1階模態(tài)分量，即有：

s1_7，重復(fù)以上步驟直至殘余分量成為一個(gè)單調(diào)函數(shù)或整個(gè)時(shí)間序列上極大值點(diǎn)數(shù)和極小值點(diǎn)數(shù)之和不大于兩個(gè)，則結(jié)束分解過程。

s2：隨機(jī)噪聲往往分布在第一個(gè)或前幾個(gè)高頻imf分量中，采用有效指標(biāo)函數(shù)eif確定各階imf分量中有效信號與噪聲的能量分界點(diǎn)，判斷含噪較多的高頻數(shù)imf分量。eif有效指標(biāo)函數(shù)定義如下：

式中，x(t)為去噪后的地震信號，n代表信號長度；imf(t)為分解出的各級模態(tài)分量；eif越小就表示imf(t)分量越近似于原始信號，eif的較大值被認(rèn)為是噪聲。

s3：依據(jù)能量分界點(diǎn)對部分高頻imf分量進(jìn)行小波閾值重構(gòu)小波系數(shù)，并保持低頻imf分量及殘余分量不變。小波閾值去噪步驟包括：

首先將能量分界點(diǎn)劃分出的含噪高頻imf分量信號轉(zhuǎn)換至小波域，求取小波系數(shù)wj：

wj＝w(imf(t))；

其中w(·)表示小波變換，j為小波分解層數(shù)；

設(shè)置一個(gè)閾值，小于該閾值的小波系數(shù)被認(rèn)為是由噪聲產(chǎn)生的，將其去除；閾值處理采用軟閾值函數(shù)

式中，sgn()為符號函數(shù)，wj為小波系數(shù)，λ為閾值，其中σ為噪聲標(biāo)準(zhǔn)差，n為信號長度；噪聲標(biāo)準(zhǔn)差median(|wj|)為小波多分辨級分解系數(shù)的中值；

最后對小波系數(shù)wj進(jìn)行逆變換，即得到小波閾值去噪后的高頻imf分量信號。

s4：將小波閾值中去噪后的高頻imf分量和不做處理的低頻imf分量以及殘余分量累加重構(gòu)，即得到去噪后的地震信號x(t)：

式中，n為ceemdan分解得到的imf分量級數(shù)；m為需要進(jìn)行小波閾值去噪的高頻imf分量數(shù)目；imfi′(t)為進(jìn)行小波閾值去噪后的高頻imf分量；imfi(t)為不需要處理的低頻imf分量；rn(t)為ceemdan分解后的殘余分量；

s5：采用改進(jìn)的滑動時(shí)窗能量比法識別地震事件并初步確定初至波走時(shí)為t0。改進(jìn)的滑動時(shí)窗能量比法的原理是根據(jù)地震記錄初至波到達(dá)時(shí)，該時(shí)刻的前后時(shí)窗能量比值a(t)達(dá)到最大值，拾取a(t)最大值對應(yīng)的時(shí)刻即為初至走時(shí)。在前、后時(shí)窗能量中加入穩(wěn)定因子，增強(qiáng)其穩(wěn)定性，改進(jìn)公式可表示為：

式中，x(t)為地震記錄振幅值，α為穩(wěn)定系數(shù)；為地震道相對能量；t1為第1個(gè)時(shí)窗起點(diǎn)，t0為第1個(gè)時(shí)窗終點(diǎn)或第2個(gè)時(shí)窗起點(diǎn)，t2為第2個(gè)時(shí)窗終點(diǎn)。

s6：去噪后地震序列可截?cái)酁槎鄠€(gè)局部平穩(wěn)部分，其各部分均可表示為自回歸ar模型；當(dāng)ar模型階數(shù)固定時(shí)，采用aic信息準(zhǔn)則法aic計(jì)算ar模型的aic函數(shù)最小值，該最小值位置，可充當(dāng)兩段地震序列的分界點(diǎn)；在對應(yīng)t0時(shí)刻前后取一時(shí)窗，使用aic函數(shù)計(jì)算在該時(shí)窗范圍內(nèi)的局部極小值，即得到精確的初至波走時(shí)t。所述去噪后地震信號的aic值采用下式計(jì)算：

式中，k為有效地震信號和隨機(jī)噪聲分量的最佳分界點(diǎn)，σ1和σ2則分別表示隨機(jī)噪聲和有效地震信號分量的方差，m表示ar模型階數(shù)，n為信號長度，c為常數(shù)；

在地震波初至來臨時(shí)波形發(fā)生較大變化，隨機(jī)噪聲和地震信號具有較大的數(shù)學(xué)統(tǒng)計(jì)差異，兩者擬合度最差，可求得最小aic值；因此，計(jì)算地震波的aic極小值可用于地震初至走時(shí)的確定。

s7：根據(jù)約束準(zhǔn)則對拾取的初至波走時(shí)進(jìn)行評價(jià)，對含異常初至?xí)r間信息的地震道進(jìn)行約束處理。根據(jù)同一炮的相鄰兩個(gè)檢波點(diǎn)的初至?xí)r間不會有較大突變，相繼各道合理的拾取時(shí)間應(yīng)形成折射段或同相軸的先驗(yàn)知識，對初至走時(shí)進(jìn)行評價(jià)，然后構(gòu)造約束準(zhǔn)則完成異常道的初至拾??；計(jì)算相鄰道的初至?xí)r間差絕對值和的平均值式中m為單炮總道數(shù)，ti是每道的初至走時(shí)，對相鄰初至?xí)r間差分進(jìn)行組合分析，如果滿足條件(-1≤n≤9,1≤i≤m,1≤j≤m)，則判定該道初至?xí)r間的可信度為1，如果不滿足，則判定該道初至?xí)r間的可信度為0；確定了部分初至?xí)r間拾取可信度為1的道，對篩選后被認(rèn)為可信度為0的道，通過道間樣條插值或者線性插值進(jìn)行重新確定，即得到異常道糾正處理后的初至走時(shí)。其中0＜k≤2，針對信噪比高的地震資料，k參數(shù)值可取較大值，反之，k取較小值。

附圖3為一種含噪地震信號初至走時(shí)拾取方法實(shí)現(xiàn)示意圖，圖3(a)為原始含噪單道地震波形，其信噪比較低(9db)，初至拾取困難；圖3(b)為本申請實(shí)施例中ceemdan聯(lián)合小波閾值去噪后地震波形，信噪比提升明顯(18db)，大部分隨機(jī)噪聲被去除，人工可準(zhǔn)確拾取地震初至?xí)r刻為492ms；圖3(c)為本申請實(shí)施例中改進(jìn)滑動時(shí)窗能量比法大致確定的地震波初至?xí)r間為486ms；圖3(d)為本申請實(shí)施例中aic信息準(zhǔn)則在486ms前后取一時(shí)窗(時(shí)窗長度為200ms)計(jì)算時(shí)窗局部范圍內(nèi)的最小aic值對應(yīng)走時(shí)為492ms，與人工拾取一致。

附圖4為合成地震道集數(shù)據(jù)去噪效果對比圖，圖4(a)為純凈的單炮記錄(基于ricker子波人工合成一組單同相軸地震數(shù)據(jù)，采樣點(diǎn)數(shù)為1000，道數(shù)為40，波速為3000m/s，道間距為50m，采樣頻率為1000hz)，圖4(b)為加入0db高斯白噪聲的含噪記錄，圖4(c)為emd去噪結(jié)果，圖4(d)為本申請實(shí)施例中ceemdan聯(lián)合小波閾值去噪結(jié)果。從圖4(c)中可以發(fā)現(xiàn)emd去噪后的記錄可以有效地將淹沒在噪聲中的同相軸恢復(fù)出來，具有一定去噪能力，但一些較小的突變?nèi)源嬖冢B(tài)混疊導(dǎo)致emd去噪過程不徹底；從圖4(d)可以看出本申請實(shí)施例的去噪效果要優(yōu)于emd方法，有效地平滑了emd方法中的高頻毛刺，使得去噪后的結(jié)果與純凈信號更為接近。

附圖5為合成單同相軸地震信號采用不同去噪方法處理后各道snr分布對比圖，由圖可知，采用本申請實(shí)施例的方法去噪后信噪比提升最為明顯。表1為各方法去噪后整體道集信噪比snr及均方根誤差rmse的數(shù)值對比，較其余方法snr整體提升明顯，重構(gòu)去噪后rmse最小。由以上圖表可知本申請實(shí)施例中聯(lián)合去噪方法可有效去除地震信號中的隨機(jī)噪聲干擾，提升地震信號信噪比，降低去噪信號的重構(gòu)誤差。

表1合成地震數(shù)據(jù)各方法去噪效果對比

其中，信噪比(snr)和均方根誤差(rmse)采用如下公式計(jì)算：

式中，si為原始地震信號，xi為去噪后地震信號，n代表信號長度。

附圖6是本申請一實(shí)施例的合成地震道集數(shù)據(jù)中各道去噪后初至拾取效果對比圖；表2為對應(yīng)拾取均方根誤差rmse及其在相應(yīng)絕對拾取誤差范圍內(nèi)的占比情況。由圖6和表2數(shù)據(jù)可知，本申請實(shí)施例中ceemdan聯(lián)合小波閾值法去除地震隨機(jī)噪聲后再采用改進(jìn)滑動時(shí)窗能量比聯(lián)合aic信息準(zhǔn)則對去噪后地震數(shù)據(jù)進(jìn)行初至拾取，其拾取誤差較小，各道初至?xí)r間分布與純凈信號經(jīng)人工拾取相當(dāng)，拾取誤差明顯低于未去噪直接采用其他初至拾取組合模式。本申請實(shí)施例中聯(lián)合去噪方法與本申請實(shí)施例中聯(lián)合初至拾取方法的組合，在50ms內(nèi)的拾取誤差范圍內(nèi)，全部準(zhǔn)確拾?。辉?0ms拾取誤差范圍內(nèi)，達(dá)95％準(zhǔn)確率；在5ms拾取誤差范圍內(nèi)，達(dá)87％準(zhǔn)確率，明顯優(yōu)于其它組合模式，驗(yàn)證了本發(fā)明方法的有效性。

表2對應(yīng)拾取誤差參數(shù)對比

附圖7是本申請實(shí)施例中四組實(shí)際單炮地震記錄，其中(a)和(b)為炸藥震源地震波形，(c)和(d)為可控源地震波形。圖中m為炮號，n為道數(shù)，dt為采樣時(shí)間間隔，dx為道間距。

附圖8為對應(yīng)四組單炮記錄采用不同初至拾取方法的拾取結(jié)果分布圖。由圖8可知本申請實(shí)施例中聯(lián)合去噪方法與本申請實(shí)施例中聯(lián)合初至拾取方法組合的拾取結(jié)果與人工拾取相當(dāng)，相比其他去噪和拾取方法組合具有較低的拾取誤差；在(c)和(d)圖中，由于無信號道或干擾道等異常情況的出現(xiàn)，導(dǎo)致采用短時(shí)窗平均/長時(shí)窗平均初至拾取方法(lta/sta)時(shí)相應(yīng)道出現(xiàn)較大拾取誤差的情況。

在本申請實(shí)施例中，根據(jù)同一炮的相鄰兩個(gè)檢波點(diǎn)的初至?xí)r間不會有較大突變，相繼各道合理的拾取時(shí)間應(yīng)形成折射段或同相軸的先驗(yàn)知識，對初至走時(shí)進(jìn)行評價(jià)，然后構(gòu)造約束準(zhǔn)則完成異常道的初至拾??；計(jì)算相鄰道的初至?xí)r間差絕對值和的平均值式中m為單炮總道數(shù)，ti是每道的初至走時(shí)，對相鄰初至?xí)r間差分進(jìn)行組合分析，如果滿足條件(-1≤n≤9,1≤i≤m,1≤j≤m)，則判定該道初至?xí)r間的可信度為1，如果不滿足，則判定該道初至?xí)r間的可信度為0；確定了部分初至?xí)r間拾取可信度為1的道，對篩選后被認(rèn)為可信度為0的道，通過道間樣條插值或者線性插值進(jìn)行重新確定，即得到異常道糾正處理后的初至走時(shí)。其中0＜k≤2，針對信噪比高的地震資料，k參數(shù)值可取較大值，反之，k取較小值。采用本申請實(shí)施例中聯(lián)合初至拾取方法可有效避免較大拾取誤差問題的發(fā)生。

圖9為對應(yīng)拾取精度roc曲線圖，roc曲線中不需要設(shè)定各方法拾取與人工拾取誤差閾值即可對算法本身效果做出評價(jià)。其中roc曲線右下方的面積被定義為auc值，其值越接近1，表示該算法準(zhǔn)確性越高。四組實(shí)際地震記錄中采用本申請實(shí)施例中聯(lián)合去噪方法與本申請實(shí)施例中聯(lián)合初至拾取方法組合的roc曲線auc平均值達(dá)0.96，明顯高于其余組合。綜合以上分析可知本申請實(shí)施例中聯(lián)合去噪方法與聯(lián)合初至拾取方法的組合模式，可取得最好的初至拾取效果。

如附圖2，本請實(shí)施例還提供一種含噪地震信號初至走時(shí)拾取裝置，包括：噪聲壓制單元u1，初至拾取單元u2，噪聲壓制單元u1與初至拾取單元u2通過電路連接。噪聲壓制單元u1中包括ceemdan分解模塊m1、imf分量分析模塊m2、小波閾值去噪模塊m3、分量重構(gòu)模塊m4，m1-m4依次按順序通過電路連接；初至拾取單元u2中包括改進(jìn)滑動時(shí)窗能量比法預(yù)拾取模塊m5、aic法再拾取模塊m6和初至評價(jià)與異常處理模塊m7，m5-m7依次按順序通過電路連接。

噪聲壓制單元u1，用于壓制原始地震信號中的隨機(jī)噪聲，是地震初至波走時(shí)準(zhǔn)確拾取的前提條件；

初至拾取單元u2，用于拾取地震波初至來臨的時(shí)間信息。

所述噪聲壓制單元u1包括：ceemdan分解模塊m1、imf分量分析模塊m2、小波閾值去噪模塊m3、分量重構(gòu)模塊m4。

ceemdan分解模塊m1，用于對原始含噪地震信號進(jìn)行逐級分解得到頻譜分布由高至低的各級imf分量以及殘余分量；

imf分量分析模塊m2，采用有效的指標(biāo)函數(shù)(eif)來確定能量分界點(diǎn)劃分分解得到的有效信號或噪聲。

小波閾值去噪模塊m3，將imf分量分析模塊m2中劃分出的含噪高頻imf分量再次進(jìn)行小波閾值去噪；

分量重構(gòu)模塊m4，用于將去噪后的高頻imf分量和低頻imf分量以及殘余量進(jìn)行累加重構(gòu)，即得到去噪后的地震信號。

所述初至拾取單元u2包括：改進(jìn)滑動時(shí)窗能量比法預(yù)拾取模塊m5、aic

法再拾取模塊m6和初至評價(jià)與異常處理模塊m7。

改進(jìn)滑動時(shí)窗能量比法預(yù)拾取模塊m5，用于識別地震事件并大致確定波至?xí)r刻，當(dāng)初至來臨時(shí)，該時(shí)刻的前后時(shí)窗能量比值a(t)將達(dá)到最大值，拾取a(t)最大值對應(yīng)的時(shí)間即為初至?xí)r間。

aic法再拾取模塊m6，用于初至?xí)r刻的精確拾取，在改進(jìn)滑動時(shí)窗能量比法預(yù)拾取模塊m5的基礎(chǔ)上，于預(yù)拾取時(shí)刻前后取一時(shí)窗，使用aic方法計(jì)算地震序列ar模型在該時(shí)窗范圍內(nèi)的局部極小值，即得到精確的初至波走時(shí)。

初至評價(jià)與異常處理模塊m7，用于處理無信號道或干擾道等異常情況出現(xiàn)時(shí)，計(jì)算機(jī)自動初至拾取方法出錯的情況。

本申請實(shí)施例一種含噪地震信號初至走時(shí)拾取裝置的各組成部分分別用于實(shí)現(xiàn)前述實(shí)施例方法的各步驟，在前述方法實(shí)施例中，已經(jīng)對各步驟進(jìn)行了詳細(xì)說明，在此不再贅述。

在一個(gè)或多個(gè)示例性的設(shè)計(jì)中，本申請實(shí)施例所描述的上述功能可以在硬件、軟件、固件或這三者的任意組合來實(shí)現(xiàn)。如果在軟件中實(shí)現(xiàn)，這些功能可以存儲與電腦可讀的媒介上，或以一個(gè)或多個(gè)指令或代碼形式傳輸于電腦可讀的媒介上。電腦可讀媒介包括電腦存儲媒介和便于使得讓電腦程序從一個(gè)地方轉(zhuǎn)移到其它地方的通信媒介。存儲媒介可以是任何通用或特殊電腦可以接入訪問的可用媒體。例如，這樣的電腦可讀媒體可以包括但不限于ram、rom、eeprom、cd-rom或其它光盤存儲、磁盤存儲或其它磁性存儲裝置，或其它任何可以用于承載或存儲以指令或數(shù)據(jù)結(jié)構(gòu)和其它可被通用或特殊電腦、或通用或特殊處理器讀取形式的程序代碼的媒介。

以上所述的具體實(shí)施例，對本申請的目的、技術(shù)方案和有益效果進(jìn)行了進(jìn)一步詳細(xì)說明，所應(yīng)理解的是，以上所述僅為本申請實(shí)施例的具體實(shí)施例而已，并不用于限定本申請的保護(hù)范圍，凡在本申請的精神和原則之內(nèi)，所做的任何修改、等同替換、改進(jìn)等，均應(yīng)包含在本申請的保護(hù)范圍之內(nèi)。