張漢羽， 劉懷山， 邢　磊??, 石旭亮，雷朝陽(yáng)

(1.中國(guó)海洋大學(xué)海底科學(xué)與探測(cè)技術(shù)教育部重點(diǎn)實(shí)驗(yàn)室， 山東 青島 266100； 2.中國(guó)科學(xué)院深?？茖W(xué)分工程研究所， 海南 三亞 572000)

?

基于偽單輸入多輸出系統(tǒng)的子波提取研究?

張漢羽1，2， 劉懷山1， 邢磊1??, 石旭亮1，雷朝陽(yáng)1

(1.中國(guó)海洋大學(xué)海底科學(xué)與探測(cè)技術(shù)教育部重點(diǎn)實(shí)驗(yàn)室， 山東 青島 266100； 2.中國(guó)科學(xué)院深?？茖W(xué)分工程研究所， 海南 三亞 572000)

摘要：?jiǎn)屋斎攵噍敵鱿到y(tǒng)提取地震子波的方法存在要求地震道之間的反射系數(shù)相同，子波多個(gè)且相異(即假設(shè)地震記錄的差異全部由子波造成，與反射系數(shù)無(wú)關(guān))，以及對(duì)噪聲敏感等缺陷。為了克服這種假設(shè)，本文利用地震數(shù)據(jù)的頻率特征，運(yùn)用二進(jìn)正交小波變換Mallat算法的思想，討論了兩種構(gòu)建偽單輸入多輸出系統(tǒng)的方法，推導(dǎo)了該系統(tǒng)下噪聲子空間提取地震子波的算法,并進(jìn)行了模型測(cè)試和實(shí)際資料驗(yàn)證工作。測(cè)試和驗(yàn)證結(jié)果顯示該方法具有較強(qiáng)的抗噪能力。

關(guān)鍵詞：PSIMO系統(tǒng)； 子空間； Mallat算法； 子波提取

引用格式：張漢羽，劉懷山，邢磊，等.基于偽單輸入多輸出系統(tǒng)的子波提取研究[J]. 中國(guó)海洋大學(xué)學(xué)報(bào)(自然科學(xué)版)， 2016,46(6): 74-83.

ZHANG Han-Yu, LIU Huai-Shan, XING Lei, et al. Based on the pseudo single-input multi-output system to extract seismic wavelet[J]. Periodical of Ocean University of China, 2016, 46(6): 74-83.

能否準(zhǔn)確地估計(jì)地震子波是地震資料波阻抗反演、反褶積處理及正演模擬成敗的關(guān)鍵因素之一[1]。目前子波提取方法主要有統(tǒng)計(jì)性方法和確定性方法[2]。統(tǒng)計(jì)性方法是指在地震子波和反射系數(shù)都未知的情況下，根據(jù)地震記錄來(lái)估計(jì)子波[2]，其不需要借助測(cè)井資料，但一般需要假設(shè)地層反射系數(shù)是高斯分布的隨機(jī)序列，子波相位是最小相位、零相位和最大相位，而實(shí)際的子波一般為混合相位[1]，或者需要運(yùn)用高階統(tǒng)計(jì)量的理論，假設(shè)反射系數(shù)是非高斯隨機(jī)序列，也要求較多的資料才能獲得較好的估計(jì)結(jié)果[2]，包括分形法[4-5]、倒雙譜法[3]、高階統(tǒng)計(jì)量FIR系統(tǒng)辨識(shí)[6-7]、同態(tài)法[8-10]、相位估計(jì)法[11-14]、全通子波匹配法[15]、EYW方程法[16-17]等。確定性方法則是通過(guò)測(cè)井資料結(jié)合井旁地震記錄來(lái)求取子波，包括貝葉斯法[18]、循環(huán)迭代法[19]、線性反演方法等。

基于單輸入多輸出系統(tǒng)(Single-Input Multi-Output System，簡(jiǎn)稱SIMO)的二階統(tǒng)計(jì)量盲信道辨識(shí)方法廣泛應(yīng)用到通訊技術(shù)等領(lǐng)域[20-21]。楊培杰等[2，22]將子波提取納入盲信號(hào)處理領(lǐng)域，提出了SIMO系統(tǒng)混合相位子波提取的算法，其不需要對(duì)反射系數(shù)分布規(guī)律和子波相位進(jìn)行假設(shè)，又具有準(zhǔn)確度高、計(jì)算量小且速度快等優(yōu)點(diǎn)，為二階統(tǒng)計(jì)量的子波提取開拓了新思路。但這種方法隱含地假設(shè)了地震道之間的反射系數(shù)相同，子波多個(gè)且相異，即地震記錄的差異全部由子波造成，與反射系數(shù)無(wú)關(guān)，且對(duì)噪聲敏感[22]。當(dāng)只有一個(gè)地震道時(shí)，也不能照搬通訊領(lǐng)域的過(guò)采樣技術(shù)將單輸入單輸出系統(tǒng)(Single-Input Single-Output System, 簡(jiǎn)稱SISO)轉(zhuǎn)化為SIMO系統(tǒng)，這使其難以在實(shí)際地震資料中應(yīng)用。

為了繼承SIMO系統(tǒng)估計(jì)子波方法，對(duì)子波相位和反射系數(shù)分布規(guī)律無(wú)假設(shè)的優(yōu)點(diǎn)，避免地震道之間反射系數(shù)相同，地震道之間的差異全部由子波造成的不合理假設(shè)，筆者從褶積模型出發(fā)，利用地震數(shù)據(jù)的頻率特征，運(yùn)用二進(jìn)正交小波變換Mallat算法的思想，討論了兩種構(gòu)建偽單輸入多輸出系統(tǒng)(Pseudo Single-Input Multi-Output System,簡(jiǎn)稱PSIMO)的方法，推導(dǎo)了PSIMO系統(tǒng)下噪聲子空間法提取子波的算法，也測(cè)試了它的抗噪能力，還應(yīng)用到實(shí)際地震資料處理中，估計(jì)出實(shí)際地震子波，并將估計(jì)子波視為已知子波，對(duì)整個(gè)剖面進(jìn)行反褶積處理。反褶積處理剖面顯示，同相軸變得更加精細(xì)、連續(xù)和平整，效果明顯。

1PSIMO系統(tǒng)的構(gòu)建

時(shí)間域地震記錄可表示為

s=w*r+v 。

(1)

式中：*為褶積運(yùn)算；s為地震記錄；w為地震子波；r為地層反射系數(shù)；v為加性高斯噪聲。將(1)式展開并離散：

(2)

式中：n為采樣點(diǎn)；L+1為子波點(diǎn)數(shù)；T為采樣間隔；J為總采樣點(diǎn)數(shù)，則離散化地震記錄的每個(gè)采樣點(diǎn)對(duì)應(yīng)的時(shí)間為t=nT。若將w視為系統(tǒng)函數(shù)，r視為激勵(lì)函數(shù)，s視為系統(tǒng)激勵(lì)響應(yīng)，則(2)式可視為一個(gè)SISO系統(tǒng)[21-22](見圖1)。

1.1 小波變換Mallat思想構(gòu)建PSIMO系統(tǒng)

A*Lo_D*Lo_R+A*Hi_D*Hi_R=A*δ=A 。

(3)

式中：A為任意層信號(hào)； Lo_D為低頻分解小波函數(shù)；Hi_D為高頻分解小波函數(shù)；Lo_R為低頻重構(gòu)小波函數(shù)；Hi_R為高頻重構(gòu)小波函數(shù)；δ為單位脈沖函數(shù)。(3)舍去抽取與插值環(huán)節(jié)的分解重構(gòu)過(guò)程，其本質(zhì)是多次褶積運(yùn)算，而褶積運(yùn)算滿足交換律，多次褶積過(guò)程可以不分先后，故對(duì)地震記錄的分解重構(gòu)等價(jià)于構(gòu)造地震子波，不影響地層反射系數(shù)(見公式(4)～(5))。

按照Mallat算法[29]二叉樹的分解重構(gòu)形式，將一個(gè)地震道s(SISO系統(tǒng)，見圖1)構(gòu)造成PSIMO(見圖2)系統(tǒng)的具體過(guò)程如下：

(1)進(jìn)行p次分解

第一次分解

u1high=s*Hi_D，

第二次分解

u2high=s*Hi_D，

……。

第p次分解

(2)將每一層分解信號(hào)重構(gòu)

第一層重構(gòu)

第二層重構(gòu)

……。

第p層低、高頻重構(gòu)

按小組進(jìn)行分組訓(xùn)練，小組同學(xué)相互交流，并對(duì)各小組學(xué)生的任務(wù)實(shí)施過(guò)程進(jìn)行全程跟蹤指導(dǎo)，并對(duì)其出現(xiàn)的問(wèn)題進(jìn)行總結(jié)和評(píng)價(jià)，指出各組學(xué)生在程序編碼過(guò)程中存在的某些共性問(wèn)題。例如，程序源文件命名不規(guī)范、程序編寫格式不規(guī)范、死循環(huán)、無(wú)注釋或注釋不清晰等問(wèn)題。

綜合上述分解和重構(gòu)過(guò)程和公式(1)，進(jìn)行p次分解時(shí)，任意層n(≤p)信號(hào)重構(gòu)表達(dá)式為

(4)

設(shè)需要構(gòu)造一個(gè)單輸入P輸出的PSIMO系統(tǒng)(見圖2)，則需獲得P個(gè)虛擬地震子道，這里稱偽子道。按照Mallat二叉樹形式，須進(jìn)行P-1次分解重構(gòu)，重構(gòu)的第i個(gè)偽子道xi可簡(jiǎn)單表示為：

(5)

(a.db4基小波等價(jià)算子時(shí)域和頻域曲線 Time domain and frequency domain curves of db4 base wavelet equivalent operator； b.中間截?cái)嗟葍r(jià)算子時(shí)域和頻域曲線 Time domain and frequency domain curves of middle truncation equivalent operator)

圖3等價(jià)算子分析

Fig.3Analysis of equivalent operator

(6)

rk=[r(k),…,r(k-N-L-q+1)]T∈(N+L+q)×1，

(7)

1.2 高信噪比小波層位構(gòu)建PSIMO系統(tǒng)

在1.1的基礎(chǔ)上，從信噪比的角度考慮：(1)地震信號(hào)和噪聲在小波層位上分布差異較大，如圖4所示，在主頻40Hz，延續(xù)時(shí)間2s的地震記錄上，加入高斯噪聲，使初始信噪比分別為50、30、20、10、0db(見圖例)，選用db1和db4基小波進(jìn)行7次分解，重構(gòu)每一層位的信號(hào)，并求得重構(gòu)信號(hào)的信噪比分布情況。從圖4a、b上可以看到，不同層位上的信號(hào)和噪聲相差較大，層位越低，信噪比越高，故可選擇高信噪比層位重構(gòu)信號(hào)；(2)SIMO系統(tǒng)盲辨識(shí)算法要求系統(tǒng)函數(shù)之間互素[25]，故在構(gòu)建PSIMO系統(tǒng)時(shí)，也須考慮構(gòu)造子波(系統(tǒng)函數(shù))的差異性。由于考慮到基小波種類多樣且互素，濾波強(qiáng)度也不一，如圖4a和b所示，db1基小波4～8層位的信噪比比原來(lái)高出2～10db，db4基小波7～8層位的信噪比比原來(lái)高出6～7db，故可選取不一樣的基小波分解信號(hào)，并只重構(gòu)回高信噪比層位構(gòu)建PSIMO系統(tǒng)，這樣既能保證構(gòu)造子波的差異性，又能提高PSIMO系統(tǒng)的信噪比。

這種構(gòu)造的具體推導(dǎo)與上文1.1相似，最終得到的表達(dá)式與(6)～(7)一樣，僅僅選用了多種基小波函數(shù)，此時(shí)等價(jià)算子可表示為

(8)

2地震子波提取

(9)

則(7)式可以寫成

(10)

(11)

則(10)式可寫成

(12)

(13)

(14)

將方法一(1.1構(gòu)造方法)和方法二(1.2構(gòu)造方法)估計(jì)的偽子波，轉(zhuǎn)化估計(jì)子波的表達(dá)式為

(15)

(16)

3模型驗(yàn)證

采用混合高斯模型，生成2000ms地層反射系數(shù)時(shí)間序列(見圖5a，僅顯示1 000ms)，采用30 Hz混合相位雷克子波(見圖6 標(biāo)號(hào)①)，兩者褶積產(chǎn)生一道地震記錄(見圖5b)，再加入一定量的白噪聲，生成含噪聲的地震記錄。

(a.地層反射系數(shù) Reflection coefficientb.合成地震記錄 Synthetic seismogram)

圖5仿真模型的建立

Fig.5The establishment of simulation model

圖6和圖7分別是方法一(1.1構(gòu)造方法)和方法二(1.2構(gòu)造方法)提取子波的結(jié)果。其中，圖6-7標(biāo)號(hào)①為時(shí)域理論子波，圖6標(biāo)號(hào)②③④⑤⑥⑦⑧⑨和圖7標(biāo)號(hào)②③④⑤分別是按照方法一和方法二構(gòu)造的理論偽子波和提取偽子波，每一個(gè)標(biāo)號(hào)對(duì)應(yīng)著兩道數(shù)據(jù)，第一道為理論偽子波，第二道為提取偽子波；圖6標(biāo)號(hào)⑩和圖7標(biāo)號(hào)⑥⑦⑧⑨分別為方法一和方法二的理論子波和提取子波，其第一道為理論子波，第二道為提取子波(見圖6a紅色曲線)。

這里須要說(shuō)明一下的是，只輸入一道地震記錄，理論子波只有一個(gè)，而這里得到4個(gè)估計(jì)子波，理論上這4個(gè)估計(jì)子波沒有任何差異，都等于理論子波，但由于偽子道的噪聲水平不一樣，可能導(dǎo)致4個(gè)子波存在差異，一般認(rèn)為與其它子波相似性較好的估計(jì)子波是準(zhǔn)確的。所以，在使用估計(jì)子波時(shí)，既可以挑選一個(gè)與其它相似性較好的子波，也可以剔除估計(jì)效果差的子波后，求取平均。圖6和圖7的a、b、c分別為信噪比為50、10和0db時(shí)子波提取的結(jié)果，圖下部分是上部分一一對(duì)應(yīng)的振幅譜。

可以看出，在信噪比為50和10db時(shí)，2種方法都能獲得較好的結(jié)果，當(dāng)信噪比為0db時(shí)，方法一結(jié)果與理論子波差別較大，而方法二子波時(shí)域形態(tài)和振幅譜都與理論子波符合較好(見圖7c紅藍(lán)線之間)，說(shuō)明本文提出的第二種方法抗噪能力較強(qiáng)，且優(yōu)于方法一。

低信噪比時(shí)，方法一提取的偽子波高頻部分會(huì)有較大缺失(見圖6c紅線和藍(lán)線之間)，這是由于重構(gòu)信號(hào)的層位越低，相對(duì)噪聲含量越大，信噪比越低(見圖4)，子波提取效果也就越差。第二種方法在選取高信噪比小波層位構(gòu)造偽子波時(shí)，既要注意子波之間的差異性，又要權(quán)衡偽子道的信噪比，見圖7b紅色曲線所示的偽子波，其與其它偽子波保持較大的差異，但主頻率比真實(shí)地震信號(hào)高很多，所以，該偽子道信噪比較低，提取效果也較差。但這種方法的性能取決于偽子波兩兩之間的最大差異，理論上所有提取的子波都是相同的。當(dāng)有個(gè)別壞道存在時(shí)，也不會(huì)影響整個(gè)子波提取的精度。

4實(shí)際地震資料應(yīng)用

圖8和9為某海域的疊后地震資料(見圖10a，僅部分顯示)采用方法一(1.1構(gòu)造方法)和方法二(1.2構(gòu)造方法)和井震聯(lián)合提取子波的結(jié)果。圖8標(biāo)號(hào)①②③④⑤⑥⑦⑧是采用db1基小波，按照方法一進(jìn)行7次分解，重構(gòu)每一層位信號(hào)(8道)組建PSIMO系統(tǒng)后，估計(jì)的偽子波結(jié)果。圖9標(biāo)號(hào)①②③④是分別采用db1、db2、db3、db4基小波，按照方法二進(jìn)行3次分解，重構(gòu)第3、4、2、4層位信號(hào)(4道)組建成PSIMO系統(tǒng)后，估計(jì)的偽子波結(jié)果。圖8標(biāo)號(hào)⑨和圖9標(biāo)號(hào)⑤⑥⑦⑧(藍(lán)色)分別是方法一和方法二提取的地震子波。圖8標(biāo)號(hào)⑩和圖9標(biāo)號(hào)⑨(紅色)為井約束條件下，確定性方法提取的子波，圖下部分為上部分一一對(duì)應(yīng)的振幅譜。

通過(guò)圖8和9對(duì)比可知，這兩種方法與井約束確定性方法，提取子波的時(shí)域波形和振幅譜都比較一致，但高頻部分，前者沒有后者豐富。圖9是將一個(gè)地震道，構(gòu)造成一個(gè)偽單輸入四輸出系統(tǒng)，而每一個(gè)偽子道都含有一個(gè)互素的偽子波，故能估計(jì)出4個(gè)偽子波和子波，但在理論上，這4個(gè)子波是一致的，都等于實(shí)際子波。圖10是方法二和井約束確定性方法提取的多個(gè)子波(圖9標(biāo)號(hào)⑤⑥⑦⑧⑨)與標(biāo)準(zhǔn)子波(令圖9標(biāo)號(hào)⑤為標(biāo)準(zhǔn)子波)相減，得到的估計(jì)子波殘差。圖10標(biāo)號(hào)①為標(biāo)準(zhǔn)子波，標(biāo)號(hào)②③④⑤⑥分別對(duì)應(yīng)圖9標(biāo)號(hào)⑤⑥⑦⑧⑨的子波殘差?？梢钥吹?，方法二提取出的4個(gè)子波，殘差都很小，標(biāo)號(hào)②④⑤幾乎為0，說(shuō)明4個(gè)提取子波幾乎無(wú)差異，與理論相符合。相對(duì)地，標(biāo)號(hào)③和⑥的殘差稍大，③是由于構(gòu)造偽子波時(shí)，重構(gòu)層位較高(db3基小波，層位2)，導(dǎo)致主頻遠(yuǎn)高于地震信號(hào)(見圖9標(biāo)號(hào)③)，信噪比相對(duì)較低，導(dǎo)致估計(jì)結(jié)果誤差較大；⑥是由于確定性子波估計(jì)方法，加入了井資料，而測(cè)井資料的分辨率高，頻帶寬，故估計(jì)子波也更加準(zhǔn)確，且分辨率高。

若提取的子波與實(shí)際情況符合，則對(duì)地震資料進(jìn)行反褶積、反演等處理，可提高地震資料的分辨能力。圖10b、c、d是分別將方法一(見圖8標(biāo)號(hào)⑨)、方法二(見圖9標(biāo)號(hào)⑧)和井震聯(lián)合法(見圖9標(biāo)號(hào)⑨)提取的子波，作為已知子波輸入，按道的順序?qū)υ嫉卣饠?shù)據(jù)(見圖11a)進(jìn)行頻率域除法反褶積后的剖面。通過(guò)圖像對(duì)比可以看到，這三種方法提取的子波反褶積后的結(jié)果，都提高了地震剖面的分辨能力：(1)使剖面的同相軸更加精細(xì)，可放大高清剖面進(jìn)行比較；(2)在一定程度上，使不連續(xù)同相軸變得更加連續(xù)。如圖11a所示，一條貫穿CDP標(biāo)號(hào)1～150的不連續(xù)同相軸(見深橙色曲線)，在藍(lán)色圓圈處存在明顯兩處中斷，但在圖11b上卻只存在一處中斷，圖c和d上不存在中斷(見放大藍(lán)色圓圈處)；(3)使不清晰、不平整或者地層太薄無(wú)法分辨的同相軸，變得清晰、平整、可分辨。圖12a、b、c、d 是圖11a、b、c、d鮮橙色矩形框內(nèi)部分剖面的放大顯示?？梢钥吹剑N方法的提取子波對(duì)剖面反褶積處理后，使原始資料某些區(qū)域不夠平整的同相軸，變得更加清晰平整(見黑色虛線標(biāo)出處)，甚至能分辨出更精細(xì)的地層(見褐色矩形框內(nèi))。

(a.地震剖面 Seismic profile； b.方法一提取子波反褶積剖面 Deconvolution profile using extracted wave of method 1； c.方法二提取子波反褶積剖面 Deconvolution profile using extracted wave of method 2； d.井震聯(lián)合法提取子波反褶積剖面Deconvolution profile using extracted wave of well seismic joint method)

圖11原始剖面和不同方法提取子波反褶積剖面對(duì)比

Fig.11Comparison between seismic section and deconvolution sections via different wavelets extracted by various methods

5結(jié)語(yǔ)

針對(duì)SIMO系統(tǒng)子波提取方法假設(shè)條件難滿足問(wèn)題，筆者利用地震信號(hào)的頻域性質(zhì)，運(yùn)用二進(jìn)正交小波分解重構(gòu)理論和Mallat算法思想，初步探討了兩種構(gòu)建PSIMO系統(tǒng)的方法，也推導(dǎo)了相應(yīng)的子波提取算法，一定程度上補(bǔ)充了SIMO系統(tǒng)提取子波的理論，也為SISO系統(tǒng)向SIMO系統(tǒng)轉(zhuǎn)化提供了方案。這兩種構(gòu)建PSIMO系統(tǒng)方法也各有一定的優(yōu)勢(shì)：第一種利用了Mallat算法能將信號(hào)分解且重構(gòu)信號(hào)的和等于初始信號(hào)的性質(zhì)，避免了反褶積運(yùn)算；第二種利用信號(hào)和噪聲在小波層位上的分布差異，構(gòu)建高信噪比的PSIMO系統(tǒng)，使提取子波得更加準(zhǔn)確。本文提出的先構(gòu)建PSIMO系統(tǒng)，再進(jìn)行子波估計(jì)算法，既不需要假設(shè)地震子波相位和反射系數(shù)分布規(guī)律，也不需要假設(shè)地震道之間反射系數(shù)的相同性，真正實(shí)現(xiàn)了無(wú)任何假設(shè)的二階累積量混合相位子波盲提取。模型試驗(yàn)和實(shí)際應(yīng)用證明了方法是有效的，且精度較高，抗噪能力強(qiáng)，具有較大應(yīng)用價(jià)值。

參考文獻(xiàn)：

[1]戴永壽, 王俊嶺, 王偉偉，等. 基于高階累積量ARMA模型線性非線性結(jié)合的地震子波提取方法研究 [J]. 地球物理學(xué)報(bào), 2008, 51(6): 1851-1859.

DAI Yong-Shou, WANG Jun-Ling, WANG Wei-Wei, et al. Seismic wavelet extraction via cumulant-based ARMA model approach with linear and nonlinear combination [J]. Chinese J Geophys (in Chinese), 2008, 51(6): 1851-1859.

[2]楊培杰, 潘勇, 穆星，等. 子空間法單輸入多輸出系統(tǒng)混合相位地震子波提取 [J]. 中國(guó)石油大學(xué)學(xué)報(bào)(自然版), 2010, 34(1): 41- 45.

YANG Pei-Jie, PAN Yong, MU Xing, et al. Mix-phase seismic wavelet extraction of SIMO system by subspace method [J]. Journal of China University of Petroleum, 2010, 34(1): 41-45.

[3]李軍紅, 張建中, 黃忠來(lái)，等. 基于倒雙譜的地震子波估計(jì)方法 [J]. 中國(guó)海洋大學(xué)學(xué)報(bào)(自然版), 2014, 44(11): 67-73.

LI Jun-Hong, ZHANG Jian-Zhong, HUANG Zhong-Lai, et al. A method for seismic wavelet estimation based on the bicepstrum of seismic data [J]. Periodical of Ocean University of China, 2014, 44(11): 67-73.

[4]樂(lè)友喜, 王才經(jīng). 任意相位子波的分形反褶積方法 [J]. 石油地球物理勘探, 1996, 31(6): 826-834.

Le Y X, Wang C J. Fractal deconvolution of any-phase wavelets [J]. Oil Geophysical Prospecting, 1996, 31(6): 826-834.

[5]趙秋亮, 李錄明, 羅省賢. 基于分形方法的地震子波提取及應(yīng)用 [J]. 石油物探, 2005, 44(1): 7-11.

ZHANG Qiu-Liang, LI Lu-Ming, LUO Sheng-Xian. Seismic wavelet extraction and application based on fractal technique [J]. Geophysical Prospecting for Petroleum, 2005, 44(1): 7-11.

[6]張賢達(dá), 保錚. 通訊信號(hào)處理 [M]. 北京：國(guó)防工業(yè)出版社，2000.

ZHANG Xian-Da, BAO Zheng. Basis for communication precess [M]. BeiJing: National Defence of Industry Press, 2000.

[7]馬建倉(cāng). 盲信號(hào)處理 [M]. 北京：國(guó)防工業(yè)出版社，2006.

MA Jian-Cang, Blind Signal Processing [M]. Beijing: National defence of Industry Press, 2006.

[8]Zhang X D，Zhang Y S. FIR system identification using high order cumulants alone. IEEE Trans on SP, 1994, 42(10): 2854-2858.

[9]Ulrych T J， Velis D R， Sacchi M D. Wavelet estimation revisited[J]. The Leading Edge, 1995, 14(11): 1139-1143.

[10]李國(guó)發(fā), 彭蘇萍, 高日勝，等. 復(fù)賽譜域提取混合相位子波的方法[J]. 天然氣工業(yè), 2005, 25(1): 85-87.

LI Guo-Fang, PENG Su-Ping, GAO Ri-Sheng, et al. Method of mixed phase wavelet extraction in cepstrum domain [J]. Natural Gas Industry, 2005, 25(1): 85-87.

[11]Petropulu A P et al. Phase reconstruvtion from bispectrum slices [J]. IEEE Trans on SP, 1998, 46(2): 527-530.

[12]Pan R L et al. Phase reconstruction in the trispectrum domain [J]. IEEETrans on Acoustics, Speech, andSP, 1987, 36(6): 895-897.

[13]李大衛(wèi), 尹成, 熊曉軍，等. 高階譜混合方法地震子波估計(jì)及處理 [J]. 地球物理學(xué)進(jìn)展, 2005, 20(1): 29-33.

LI Da-Wei, YI Cheng, XIONG Xiao-Jun, et al. A hybrid high order spectrum method for seismic wavelet estimation and processing [J]. Progress in Geophysics, 2005, 20(1): 29-33.

[14]謝桂生, 石玉梅，魏野. 雙譜地震子波估計(jì) [J]. 西南石油學(xué)院學(xué)報(bào), 2000, 22(3): 25-28.

XIE Gui-Sheng, SHI Yu-Mei, WEI Ye. Wavelet estimation based on bispectrum [J]. Journal of Southwest Petroleum Institute, 2000, 22(3): 25-28.

[15]Misra S, Sacchi M. Wavelet estimation by non-linear optimization of all pass operators[J]. CSEG RECORDER. 2006, 38-42.

[16]Porsani M J, Vrsin B. Mixed-phase deconvolution [J]. Geophysics, 1998, 63(2): 637-647.

[17]Ursin B, Porsani M J. Estimation of an optimal mixed-phase inverse filter [J]. Geophysical Prospecting, 2000, 48(4): 663-676.

[18]Buland Omre H. Bayesian wavelet estimation from seismic and well data [J]. Geophysics, 2003, 68(6): 2000-2009.

[19]張廣智, 劉洪, 印興耀. 井旁道地震子波精細(xì)提取方法 [J]. 石油地球物理勘探, 2005, 40(2): 158-162.

ZHANG Guang-Zhi, LIU Hong, YIN Xing-Yao. Method for fine picking up seismic wavelet at uphole trace [J]. Oil Geophysical Prospecting, 2005, 40(2): 158-162.

[20]Lang Tong. Blind identification and equalization based on second-order statistics: a time domain approach [J]. IEEE Transacticins on Information Theory, 1994, 40(2)： 340-349.

[21]蔣靜. 基于子空間的二階統(tǒng)計(jì)量盲信道辨識(shí)算法研究[D]. 鄭州：鄭州大學(xué), 2010.

JIANG Jing. Blind Channel Identification Based on Subspace Method of Second-order Statistics[D]. Zhengzhou：Zhengzhou University, 2010.

[22]楊培杰. 地震子波盲提取與非線性反演[D]. 北京：中國(guó)石油大學(xué), 2008.

YANG PeiJie. Seismic Wavelet Blind Extraction and Non-Liner Inversion[D]. Beijing： China University of Petroleum(East China), 2008.

[23]邱娜. 地震子波分解與重構(gòu)技術(shù)研究[D]. 青島：中國(guó)海洋大學(xué), 2012.

QIU Na. Research on Seismic wavelet Decomposition and Reconstruction Technology[D]. Qingdao： Ocean University of China, 2012.

[24]Eric Moulines, Pierre Duhamel. Subspace Methods for the Blind Identification of Multichannel FIR Filters [J]. IEEE Transactions on Signal Processing, 1995, 43(2)： 516-525.

[25]何子蘇, 夏威. 現(xiàn)代數(shù)字處理及其應(yīng)用 [M]. 北京：清華大學(xué)出版社, 2004.

HE Zi-Su, XIA Wei. Modern digital processing and its application [M]. Beijing: Tsinghua University Press, 2004.

[26]Eric Moulines, Pierre Duhamel. Subspace methods for the blind identification of multichannel FIR filters [J]. IEEE Transactions on Signal Processing, 1995, 43(2)：516-525.

[27]Xia W, He Z H, Liu B Y, et al. Adaptive multichannel blind identification using manifold optimization [J]. Signal Processing, 2008(88): 1595-1605.

[28]周偉. 基于MATLAB的小波分析應(yīng)用(第2版) [M]. 西安：西安電子科技大學(xué)出版社出版, 2010.

ZHOU Wei. Wavelet analysis applications based on MATLAB(second edition)[M]. XiAn: XiAn Electronic Science & Technology University Press, 2010.

[29]葛哲學(xué), 沙威. 小波分析理論與MATLAB R2007實(shí)現(xiàn) [M]. 北京: 電子工業(yè)出版社, 2007, 10: 76-85.

GE Zhe-Xue, SHA Wei. Wavelet analysis theory and implementation of R2007 MATLAB[M]. BeiJing: Electronic Industry Press, 2007, 10: 76-85.

責(zé)任編輯徐環(huán)

Based on the Pseudo Single-Input Multi-Output System to Extract Seismic Wavelet

ZHANG Han-Yu1,2, LIU Huai-Shan1, XING Lei1, SHI Xu-Liang1, LEI Chao-Yang1

(1.Key Lab of Submarine Geosciences and Prospecting Techniques, Ministry of Education, Ocean University of China, Qingdao 266100, China; 2.Institute of Deep-sea Science and Engineering, Chinese Academy of Sciences, Sanya 572000, China)

Abstract:The method of extracting seismic wavelet based on the Single-input Multi-output (SIMO) System, requires reflection coefficients to be the same as the others and wavelets to be different from the others. In other words, it unreasonably assumes that all the difference of seismic record is caused by multi-wavelet and has nothing to do with reflection coefficients. In order to avoid these assumptions, Aiming at the character of seismic signal and using the thinking about Mallat algorithm of the Discrete Wavelet Transform, this paper makes a preliminary discussion about two methods of building Pseudo Single-Input Multi-Output System (PSIMO), and then deduces the algorithm of extracting seismic wavelet with noise subspace method in PSIMO system. The results of simulation and application demonstrate that those methods are effective and characterized by strong stability to resist noise.

Key words:PSIMO； subspace； mallat algorithm； wavelet extraction

基金項(xiàng)目：? 國(guó)家自然科學(xué)基金項(xiàng)目(41176077 、41230318和41304096)；教育部博士點(diǎn)基金項(xiàng)目(20130132120014)；國(guó)家高技術(shù)研究與發(fā)展計(jì)劃項(xiàng)目(2013AA092501)資助

收稿日期：2015-02-03；

修訂日期：2015-04-05

作者簡(jiǎn)介：張漢羽(1989-)，男，碩士生，主要從事海洋地震勘探研究。E-mail: 824284715@qq.com ??通訊作者：E-mail: xingleiouc@ouc.edu.cn

中圖法分類號(hào)：P315.3

文獻(xiàn)標(biāo)志碼：A

文章編號(hào)：1672-5174(2016)06-074-10

DOI：10.16441/j.cnki.hdxb.20150029

Supported by The National Natural Science Foundation of China(41176077、41230318、41304096)、Ph.D. Programs Foundation of Ministry of Education of China(2013132120014) and The National High Technology Research and Development Program of China(2013AA092501).