多通道瞬變電磁m序列全時正演模擬與反演

齊彥福, 殷長春*， 王若, 蔡晶

1 吉林大學(xué)地球探測科學(xué)與技術(shù)學(xué)院， 長春 130021 2 中國科學(xué)院地質(zhì)與地球物理研究所， 北京 100029

?

多通道瞬變電磁m序列全時正演模擬與反演

齊彥福1, 殷長春1*， 王若2, 蔡晶1

1 吉林大學(xué)地球探測科學(xué)與技術(shù)學(xué)院， 長春 130021 2 中國科學(xué)院地質(zhì)與地球物理研究所， 北京 100029

傳統(tǒng)瞬變電磁方法主要用于金屬礦勘查，無法滿足油氣資源高阻目標(biāo)體的勘探需要.多通道瞬變電磁(MTEM)系統(tǒng)的出現(xiàn)解決了這一問題.該方法采用偽隨機序列發(fā)射波形和擬地震觀測方式，測量同線電場分量，記錄全時發(fā)射電流及多道觀測數(shù)據(jù)，實現(xiàn)對高阻薄層的高精度探測.鑒于國內(nèi)對此方法的研究還處于理論探索階段，尚未進行相應(yīng)的仿真模擬和數(shù)據(jù)處理工作，本文針對m序列發(fā)射波形多通道瞬變電磁法的全時正演模擬和反演解釋進行研究，為國內(nèi)正在進行的MTEM儀器系統(tǒng)研發(fā)及數(shù)據(jù)解釋提供理論指導(dǎo).我們利用方波響應(yīng)移位疊加和電流導(dǎo)數(shù)與階躍響應(yīng)褶積兩種方法實現(xiàn)理論m序列和實際發(fā)射波形的全時正演模擬；再通過相關(guān)辨識技術(shù)，削弱噪聲影響，計算脈沖響應(yīng)；最后對積分得到的階躍響應(yīng)進行共中心點道集數(shù)據(jù)聯(lián)合反演，獲取地下電性分布信息.

多通道瞬變電磁； m序列； 全時正演； 反演

The MTEM system adopts quasi-seismic observation, measuring in-line electric field component of pseudo-random binary sequence (PRBS) waveforms. Full-wave transmitting current and multi-channel observation data are recorded, so that the resistive thin layers can be detected. We perform full-wave forward modeling for a theoretical m-sequence and actual transmitter waveforms based on superposition of square wave responses and convolution of step responses with derivative of transmitting current, respectively. Taking advantages of the correlation identification technique, the impulse response is calculated while the noise is suppressed. Finally, we invert the CMP data by the Occam′s algorithm to obtain underground electrical information.

Firstly, Fourier transform is used to analyze the spectral components of square waveforms, 5-order 2n-sequence and 5-order m-sequence. It shows that the m-sequence has a larger frequency bandwidth with an equal interval between spectral lines than another two waveforms. Secondly, in order to verify the full-wave forward modelling results, a homogenous half-space of 30 Ωm is set up with a suit of 4-order m-sequence transmitting current. The results obtained by superposition of square wave responses, convolution of step responses with derivate of transmitting current and convolution of impulse response with transmitting current are compared. It shows that the first two results are consistent. We further add 30 dB Gaussian noise to the responses. The relative error after correlation between modeling response and transmitting waveforms is reduced to 1/4 of the original one. The step response obtained by integrating the impulse response is very close to the theoretical value with relative error of just 0.08%. Finally, we design two three-layer models to demonstrate the capability of the MTEM method to resolve shallow and deep hydrocarbon reservoirs. In the first model, a reservoir of 500 Ωm and 50 m thickness is set at the depth of 300 m in a homogenous half-space of 50 Ωm. The survey is done respectively at offsets of 900, 1000, 1100, 1200 and 1300 m. For the second model, we assume the hydrocarbon reservoir is moved to 2000 m depth, the survey is done at 5000, 5200, 5500, 5700 and 6000 m offsets, respectively. The inversion results for the two models by the Occam′s method show that the MTEM method can well resolve both shallow and deep hydrocarbon reservoirs. The CMP data have better resolution than one-offset data.

The MTEM method using m-sequence as transmitting waveforms, measuring multi-channel in-line electric fields, turns out to have a good resolution to hydrocarbon reservoirs. Spectral analysis shows that the m-sequence has a wide frequency band with an equal interval between spectrum lines. The techniques based on superposition of square-wave responses and convolution of step responses with derivate of transmitting current are proved to be very effective for modeling theoretical m-sequence and actual full-wave responses. Using correlation identification technique, the influence from transmitting waveforms is removed in the calculation of impulse response while the noise is suppressed. The Occam′s inversion of CMP data obtained by integration of impulse response can well resolve underground resistive targets like reservoirs.Keywords Multi-transient electromagnetic (MTEM); m-sequence; Full-wave forward modeling; inversion

1 引言

近年來，電磁方法在資源勘探中獲得越來越廣泛的應(yīng)用.瞬變電磁法作為一種重要的電磁勘探手段越來越受到重視.然而，由于電磁干擾的影響以及目標(biāo)體埋深的不斷加大，對瞬變電磁法勘探能力的要求不斷提高，新方法、新技術(shù)的研究成為解決這一問題的有效途徑(薛國強等，2013).

考慮到發(fā)射波形對瞬變電磁法勘探能力的影響，前人對傳統(tǒng)的方波、梯形波、半正弦波和三角波等發(fā)射波形進行了深入研究(Liu，1998；陳曙東等，2012；關(guān)珊珊等，2012)，以期在原有裝置和觀測系統(tǒng)基礎(chǔ)上，通過調(diào)整發(fā)射波形及相關(guān)參數(shù)達到改善勘探效果的目的.為獲得更加豐富的頻率信息，何繼善(1998)提出了偽隨機脈沖電磁觀測方法和技術(shù).在儀器設(shè)備、裝置型式和數(shù)據(jù)處理手段上前人也做了許多改進.陳曉東等(2012)研制出適合瞬變電磁法的高溫超導(dǎo)磁強計，利用其直接測量磁場分量，與傳統(tǒng)感應(yīng)線圈傳感器相比有效地提高了探測深度.李貅等(2012)將合成孔徑應(yīng)用到瞬變電磁數(shù)據(jù)處理當(dāng)中，利用相關(guān)疊加技術(shù)，實現(xiàn)多孔徑數(shù)據(jù)合成，將單點處理方式發(fā)展成為逐點推移多次覆蓋的處理方法，大大提高了瞬變電磁法的分辨率.薛國強等(2013)提出短偏移距瞬變電磁技術(shù)(SOTEM).隨著偏移距的縮短，接地導(dǎo)線源的場對地層的反映更加靈敏，通過增加觀測時長，加大探測深度，使該方法成為大深度、高分辨探測地下礦產(chǎn)資源的一種新技術(shù)手段.

Wright 等(2002)設(shè)計出一種全新的電磁觀測系統(tǒng)——多通道瞬變電磁系統(tǒng)(MTEM)，此方法利用偽隨機序列作為發(fā)射波形，測量同線電場分量，同時觀測多道數(shù)據(jù).Wright將其成功應(yīng)用于油氣資源勘查和油田開采狀態(tài)監(jiān)測，大大拓寬了瞬變電磁法的應(yīng)用領(lǐng)域.Ziolkowski 等(2007)利用MTEM裝置對法國西南部天然氣儲層進行測量，獲得了理想的勘探效果.Wright 和Ziolkowski(2007)分析了MTEM數(shù)據(jù)噪聲的類型及相應(yīng)壓制方法.Ziolkowski 等(2010)再次對北海Harding油氣田進行了成功勘探.Ziolkowski等(2011)在挪威Peon 油氣田利用瞬變電磁方法進行淺?？碧?，并將方波和偽隨機兩種信號源進行了對比，結(jié)果顯示偽隨機信號擁有更強的分辨能力.

m序列偽隨機碼作為目前應(yīng)用最廣泛的一種偽隨機系列已被應(yīng)用于通信領(lǐng)域，如擴頻通信、衛(wèi)星通信的碼分多址、數(shù)字?jǐn)?shù)據(jù)中的加密等領(lǐng)域(林智慧等，2009).鑒于m序列具有諸多優(yōu)良的性質(zhì)，逐漸受到地球物理工作者的關(guān)注.湯井田和羅維斌(2008)對基于相關(guān)辨識的逆重復(fù)m序列偽隨機電磁法的原理進行分析，討論參數(shù)選取對壓制噪聲、提高信噪比的影響.逆重復(fù)m序列是由m序列隔位取反得到的，它們具有十分相近的性質(zhì)，對m序列的研究具有重要的參考價值.Li 等(2012)對相關(guān)識別技術(shù)在時間域譜激電中的應(yīng)用做了更加深入的研究.淳少恒等(2014)對偽隨機m序列及其在電法勘探中的應(yīng)用現(xiàn)狀做了比較全面的介紹.

與國外日趨成熟的技術(shù)相比，國內(nèi)對m序列瞬變電磁法的研究還停留在理論探索階段，并沒有進行仿真模擬和數(shù)據(jù)處理研究(淳少恒等，2014).本文在前人已有工作基礎(chǔ)上，將現(xiàn)有技術(shù)集成應(yīng)用到一個新的領(lǐng)域——多通道瞬變電磁法，并給出了針對MTEM裝置m序列發(fā)射波形正演模擬及反演解釋方法和完整技術(shù)流程，為國內(nèi)正在進行的MTEM系統(tǒng)儀器研究及其后續(xù)數(shù)據(jù)解釋工作提供理論指導(dǎo).我們首先簡要介紹m序列的產(chǎn)生及其重要性質(zhì)，然后著重闡述m序列響應(yīng)的全時正演模擬以及脈沖辨識方法，最后對積分獲得的階躍響應(yīng)進行Occam反演，獲得地下電性分布信息.本文的研究成果必將對國內(nèi)m序列多通道瞬變電磁勘查系統(tǒng)研發(fā)起到積極的促進作用.

2 MTEM裝置及測量方式

根據(jù)Ziolkowski等(2007)，MTEM裝置采用電偶極發(fā)射源，通過電極AB向地下發(fā)射偽隨機序列(PRBS)，在同線方向觀測電場ex分量(如圖1).在整個工作過程中，觀測并記錄全部時刻的發(fā)射電流和電磁響應(yīng).

圖1 MTEM裝置示意圖Fig.1 Schematic diagram of MTEM

在實際工作中，一般連續(xù)發(fā)射多個周期的偽隨機序列，以增大信號強度.對于同一套偽隨機編碼序列進行重復(fù)發(fā)射，可以實現(xiàn)觀測數(shù)據(jù)的多次疊加，從而壓制噪聲.通過移動發(fā)射源和接收機位置，可以獲得共中心點道集多道觀測數(shù)據(jù).對其同時反演處理，可以有效地提高對高阻儲集層的識別能力.

3 m序列

3.1 產(chǎn)生機理

m序列是由線性反饋移位寄存器(圖2)生成的最長線性反饋移位寄存器序列，其周期為2n-1，n為移位寄存器的階數(shù).通過調(diào)整移位寄存器起始狀態(tài)和反饋線鏈接狀態(tài)，可以實現(xiàn)對m序列進行編碼(吳明捷和杜天蒼，2002).

圖2 線性反饋移位寄存器示意圖Fig.2 Schematic diagram of linear feedback shift register

圖2為線性反饋移位寄存器示意圖，其中a為移位寄存器的狀態(tài)，C為反饋系數(shù)，C=0表示斷開，而C=1表示鏈接.碼元的新狀態(tài)由n階串接寄存器線性抽頭經(jīng)模2加法反饋獲得(淳少恒等，2014).

3.2 m序列特性

m序列具有很多優(yōu)良性質(zhì)：(1)均衡性：在一個周期內(nèi)“1”元素出現(xiàn)的次數(shù)比“0”元素出現(xiàn)的次數(shù)多1；(2)游程特性：一個周期內(nèi)共有2n-1個游程，且“0”元素的游程和“1”元素的游程各占一半，長度為k(1≤k≤n-2)的游程數(shù)占游程總數(shù)的1/2k(林智慧等， 2009)；(3)優(yōu)良的自相關(guān)特性：m序列具有尖銳的自相關(guān)函數(shù)，其表達形式為(謝雪康和楊曉蓉，2010)

(1)式中N是m序列碼長，t0是碼元寬度，k為周期數(shù)(k=0,1,2,…).由(1)式可以看出碼長越長且碼元寬度越小，m序列的自相關(guān)函數(shù)越接近脈沖函數(shù).

3.3 頻譜成分分析

在傳統(tǒng)電磁法中，一般發(fā)射方波作為激勵源.這種波形在電子元件上容易實現(xiàn)，然而其頻率成分比較單一，能量主要集中在基頻及其奇次諧波上(Ziolkowski et al.，2011)，且隨頻率增加能量迅速衰減，基波及3次諧波的功率占方波總功率的90%以上(湯井田和羅維斌，2008).為了豐富發(fā)射波形的頻譜成分，何繼善(1998)提出了偽隨機觀測方案，設(shè)計并發(fā)射2n系列偽隨機信號(何繼善等，2009).2n系列的能量主要集中在主頻2k(k=0，1，2，…，n-1)上，且主頻上的幅值基本相等，隨著n的增大頻帶變寬(湯井田和羅維斌，2008).與前兩種發(fā)射波形不同，m序列偽隨機碼具有較寬的頻帶，譜線之間是等間隔的，包含更多的頻率信息，而且其包含直流成分，這是前兩種波形所沒有的(湯井田等，2007).

將m序列作為發(fā)射電流波形，影響其勘探能力的主要參數(shù)包括m序列的階數(shù)、碼元寬度、重復(fù)發(fā)射周期數(shù)以及電流強度(湯井田和羅維斌，2008；淳少恒等，2014).碼元寬度和階數(shù)影響頻帶寬度及能量分布.碼元寬度越大，階數(shù)越低，頻帶就越窄，能量主要集中在低頻部分，勘探深度越大，信噪比也越高.然而由于高頻成分的缺失必將造成分辨率降低.碼元寬度越小，階數(shù)越高，頻帶就越寬，高頻成分增加，分辨率將得到提高.然而由于能量被均勻分配到整個頻帶范圍內(nèi)，勢必造成單個頻率成分能量降低，降低信噪比.電流強度越大，重復(fù)發(fā)射次數(shù)越多，壓制噪聲的能力就越強，信噪比越高.隨之對發(fā)射機的硬件要求就越高，工作效率降低.在實際工作當(dāng)中，應(yīng)綜合考慮信噪比、分辨率和勘探深度，根據(jù)勘探目標(biāo)合理調(diào)整m序列的相關(guān)參數(shù)，才能達到理想的勘探效果.

圖3 發(fā)射波形及頻譜特征 (a1) 方波波形；(b1) 5階2n序列；(c1) 5階m序列； (a2)方波振幅譜；(b2) 5階2n序列振幅譜；(c2) 5階m序列振幅譜.Fig.3 Transmitted waveform and spectrum (a1) Square wave; (b1) 5 series 2n-sequence; (c1) 5 series m-sequence; (a2) Amplitude spectrum for square wave; (b2) Amplitude spectrum for 5 series 2n-sequence; (c2) Amplitude spectrum for 5 series m-sequence.

4 m序列全時正演理論

4.1 理論m序列電磁響應(yīng)模擬

在數(shù)字信號處理中，方波可以利用兩個階躍波反向移位疊加合成(圖4).參照此原理，在忽略二次電流耦合的條件下，可采用階躍響應(yīng)反向移位疊加計算方波響應(yīng)(圖5).

圖4 階躍波合成方波示意圖Fig.4 Schematic diagram of square-wave obtained from superposed step-waves

圖5 階躍響應(yīng)合成方波響應(yīng)Fig.5 Square-wave response obtained from step-waves

圖5中帶有圓圈符號的灰色實線為第一個正向階躍波激發(fā)的電磁響應(yīng)，帶有三角符號的灰色實線是第二個負(fù)向階躍波激發(fā)的電磁響應(yīng)，黑色實線是兩個階躍響應(yīng)疊加合成的方波響應(yīng).

m序列實際上是由一系列方波組成的，所以m序列可以利用方波移位疊加進行合成(圖6).與階躍響應(yīng)合成方波響應(yīng)的方法類似，可采用方波響應(yīng)合成m序列響應(yīng).

圖6 方波合成m序列示意圖Fig.6 Schematic diagram of m-sequence obtained from square-waves

4.2 實際發(fā)射波形電磁響應(yīng)模擬

在實際工作中，發(fā)射電流在供電瞬間并不能迅速達到穩(wěn)定，斷電時出現(xiàn)延遲，發(fā)射波形實際上為一系列梯形波.為模擬實際發(fā)射波形的電磁響應(yīng)，可以采用褶積方法進行計算.Smith 和Lee(2004)指出：在航空系統(tǒng)中，當(dāng)t→0時半空間的脈沖響應(yīng)趨近于t-1/2，出現(xiàn)奇異性，而階躍響應(yīng)趨近一個常數(shù).針對這一問題，殷長春等(2013)利用階躍響應(yīng)與電流的時間導(dǎo)數(shù)褶積代替脈沖響應(yīng)與電流褶積的方法成功實現(xiàn)了任意波形全時響應(yīng)的計算.研究發(fā)現(xiàn)，在地面電磁法中，t→0時電場脈沖響應(yīng)同樣出現(xiàn)奇異性，所以本文采用階躍響應(yīng)與電流導(dǎo)數(shù)褶積的方法計算實際發(fā)射波形電磁響應(yīng)，即

(2)

其中I(t)是發(fā)射電流，而Bs(t)是下階躍響應(yīng).(2)式中只考慮了由于電流變化感應(yīng)出的電磁場，并沒有考慮電流自身產(chǎn)生的直流場，為計算任意發(fā)射波形全時響應(yīng)的總場，需在(2)式中加入直流項，

(3)

式中fDC是單位電流強度的直流場.

4.3 結(jié)果對比

理論上，當(dāng)發(fā)射波形供電穩(wěn)定和斷電延遲時間足夠短時，梯形波可以近似方波，利用移位疊加和褶積兩種方法計算的m序列電磁響應(yīng)具有可比性.本文設(shè)計一個30 Ωm的均勻半空間模型，電偶源采用單位長度，收發(fā)距為1000 m，電流強度為30 A，發(fā)射波形為4階m序列(圖7)，正負(fù)電平轉(zhuǎn)換時間為40.96 μs.

圖7 4階m序列電流波形Fig.7 Transmitting waveform of 4-order m-sequence

圖8 移位疊加與褶積結(jié)果對比圖Fig.8 Comparison of shift-superposition with convolution results

圖8為電場ex響應(yīng)，其中黑色實線為方波響應(yīng)移位疊加計算結(jié)果，十字星為階躍響應(yīng)與電流導(dǎo)數(shù)褶積計算結(jié)果，而灰色實線為脈沖響應(yīng)與電流波形褶積結(jié)果.可以看出，三種方法在off-time會得到一致的結(jié)果，然而在on-time，由于脈沖響應(yīng)在t→0時出現(xiàn)奇異性，使得結(jié)果出現(xiàn)明顯偏差.而前兩組結(jié)果吻合較好，均方相對誤差小于3%.

5 m序列電磁響應(yīng)反演

由于在電磁反演當(dāng)中需要多次重復(fù)正演運算，

如果直接對含有發(fā)射波形的電磁響應(yīng)進行反演，正演計算量過大，嚴(yán)重影響反演效率.因此本文利用反卷積方法，消除發(fā)射波形的影響，提高反演效率.

5.1 脈沖辨識

根據(jù)Ziolkowski等(2007)，m序列的電磁響應(yīng)可以寫為

(4)

式中I(t)是發(fā)射電流，BI(t)是脈沖響應(yīng)，rCD(t)是接收機響應(yīng)，n(t)是噪聲.通常情況下，將發(fā)射電流和接收機響應(yīng)的褶積稱為儀器系統(tǒng)響應(yīng)，可將(4)式化簡為

(5)

式中s(t)是儀器系統(tǒng)響應(yīng)，具有與發(fā)射電流相同的性質(zhì).根據(jù)Wiener-Hopf方法(Wiener, 1949)，將(5)式兩端分別與儀器系統(tǒng)響應(yīng)做互相關(guān)，可以得到

(6)

其中Rss(t)是s(t)的自相關(guān)矩陣.由于m序列具有優(yōu)良相關(guān)特性，噪聲是隨機的，所以n(t)?s(t)≈0，(6)式可化簡為

(7)

式中B′(t)=B(t)?s(t).當(dāng)Rss(t)為脈沖函數(shù)時可以得到

(8)

但實際工作中Rss(t)并不是脈沖函數(shù)，所以將(8)式寫成離散矩陣形式

(9)

(9)式可采用最小二乘方法進行求解.

圖9展示了利用m序列的相關(guān)特性壓制噪聲的效果.圖9a為電場響應(yīng)ex，其中黑色實線是理論電場響應(yīng)，灰色實線是加入高斯噪聲后電場響應(yīng)，其信噪比為30 dB.圖9b為電場ex響應(yīng)與發(fā)射波形(圖7)互相關(guān)后的結(jié)果，其中黑色實線和灰色實線分別為圖9a中理論和加入噪聲后響應(yīng)與發(fā)射波形的互相關(guān)結(jié)果，兩組結(jié)果幾乎完全重合，而均方相對誤差減少1/4多.可以看出利用m序列的相關(guān)特性可以有效地壓制噪聲干擾.

本文針對圖8中的m序列響應(yīng)進行相關(guān)辨識，獲得脈沖響應(yīng)辨識結(jié)果(圖10)，其中黑色實線為理論脈沖響應(yīng)，十字星為相關(guān)辨識結(jié)果，可以看出在早期兩組結(jié)果吻合較好；在晚期(20 ms之后)，辨識結(jié)果出現(xiàn)微弱震蕩.

5.2 積分計算階躍響應(yīng)

根據(jù)裴易峰等(2014)，采用簡單的矩形數(shù)值積分算法可以實現(xiàn)由dB/dt計算B，通過對晚期道數(shù)據(jù)進行指數(shù)擬合可以有效壓制噪聲.

(10)

B(0)是上階躍響應(yīng)的起始值，為0.針對圖10中利用相關(guān)辨識方法獲得的脈沖響應(yīng)，采用此方法計算階躍響應(yīng).

圖9 相關(guān)壓制噪聲效果Fig.9 Noise suppression via correlation

圖10 脈沖辨識結(jié)果Fig.10 Impulse response identified by deconvolution

圖11是理論階躍響應(yīng)與積分結(jié)果對比圖，其中黑色實線是理論階躍響應(yīng)，十字星是積分結(jié)果.可以發(fā)現(xiàn)兩組結(jié)果十分吻合，平均相對誤差為0.08%.積分運算有效地壓制了辨識脈沖響應(yīng)晚期震蕩產(chǎn)生的影響.

5.3 數(shù)據(jù)反演

本文對積分獲得的階躍響應(yīng)利用Occam方法(Constable et al., 1987)進行反演，此方法在目標(biāo)函數(shù)中加入光滑度約束，使得在滿足數(shù)據(jù)擬合的同時，獲得最光滑的反演模型.其目標(biāo)函數(shù)定義為

(11)

式中?是粗糙度矩陣，m是模型參數(shù)，μ是拉格朗日乘子，W是均方差矩陣，F(xiàn)(m)是正演響應(yīng).通過求目標(biāo)函數(shù)的極值，可以獲得新模型

(12)

圖11 理論階躍響應(yīng)與積分結(jié)果對比Fig.11 Comparison of theoretical step response with integrated result

其中J是雅克比矩陣.采用迭代技術(shù)，對反演模型進行更新.在每次迭代過程中搜索合適的拉格朗日乘子，調(diào)整光滑約束項和數(shù)據(jù)擬合項的相對權(quán)重，使反演模型逐漸逼近真實解，最終獲得最光滑的反演模型.此方法反演穩(wěn)定，受初始模型影響小.

6 算例

為了檢驗MTEM系統(tǒng)對高阻儲集層的勘探能力，本文設(shè)計了儲層埋深不同的兩個模型.針對淺部高阻模型，發(fā)射高階、碼元寬度較小的m序列，增加高頻發(fā)射成分，提高分辨率；對于深部高阻模型，發(fā)射低階、碼元寬度較大的m序列，增強低頻發(fā)射成分，加大勘探深度.

6.1 模型一

為檢驗MTEM裝置對淺部高阻儲層的反映能力，本文首先設(shè)計了一個三層介質(zhì)模型(圖12)，在50 Ωm的均勻半空間中存在500 Ωm的高阻儲集層，埋深為300 m.發(fā)射源長度為50 m，發(fā)射波形為8階m序列偽隨機碼，碼元寬度為1 ms(圖13).根據(jù)Ziolkowski等(2007)，當(dāng)收發(fā)距大于2倍儲集層埋深時，脈沖響應(yīng)的峰值強度和出現(xiàn)時刻會出現(xiàn)明顯的異常.基于這一結(jié)論，本文選擇了900，1000 m，1100 m，1200 m和1300 m五組收發(fā)距進行觀測.

圖14是收發(fā)距為1000 m的m序列電場響應(yīng).其中黑色實線為不含有高阻層的均勻半空間響應(yīng)，灰色實線為含有高阻儲層的電磁響應(yīng).由圖可以看

圖12 淺部高阻儲集層模型Fig.12 Model for a shallow hydrocarbon reservoir

圖13 8階m序列發(fā)射波形Fig.13 8-order m-sequence transmitting waveform

出高阻儲層對電場的影響非常大.在on-time和off-time與均勻半空間模型均有較大差別，在off-time出現(xiàn)一個較高的負(fù)向峰值.在正演數(shù)據(jù)中加入高斯噪聲，信噪比為65 dB，然后利用相關(guān)辨識技術(shù)獲得脈沖響應(yīng)，并積分計算階躍響應(yīng).圖15為淺層高阻儲層存在與否對脈沖和階躍響應(yīng)影響對比圖.圖15a是利用相關(guān)辨識技術(shù)獲得的脈沖響應(yīng)，圖15b是對圖15a中脈沖響應(yīng)進行積分計算的階躍響應(yīng).其中黑色實線為不含高阻儲層的均勻半空間響應(yīng)，而灰色實線為含有高阻儲層的響應(yīng)曲線.由圖可以看出，由于噪聲的存在使得脈沖響應(yīng)出現(xiàn)了震蕩，而積分方法有效地壓制了噪聲干擾，獲得光滑的階躍響應(yīng)曲線.與均勻半空間情況相比，高阻儲層的出現(xiàn)使得脈沖響應(yīng)出現(xiàn)較高的峰值，而且峰值時間向早期道偏移.除此之外，階躍響應(yīng)的振幅也明顯增強.

通過移動發(fā)射源位置，獲得多組共中心點道集(CMP)數(shù)據(jù)，利用Occam方法進行反演.圖16為不同收發(fā)距的反演結(jié)果，其中黑色粗實線、黑色細實線、灰色粗實線、黑色虛線和灰色細實線分別是收發(fā)距為900、1000、1100、1200 m和1300 m的反演結(jié)果，灰色虛線是對以上5個收發(fā)距數(shù)據(jù)同時反演結(jié)果.可以看出單個收發(fā)距數(shù)據(jù)的反演結(jié)果均可以反映出地下電性分布的總體趨勢，然而對高阻儲層位置和電阻率的反映不是很準(zhǔn)確.與其相比，利用5個收發(fā)距數(shù)據(jù)同時反演獲得的結(jié)果更加接近真實模型，界面位置和高阻儲層電阻率都比較準(zhǔn)確.因此，利用MTEM裝置同時測量多道數(shù)據(jù)探測淺部高阻儲層十分有效.

6.2 模型二

為檢驗MTEM裝置對深部高阻儲層的反映能力，本文再次設(shè)計了一個三層介質(zhì)模型(圖17)，在50 Ωm的均勻半空間中存在500 Ωm的高阻儲集層，埋深為2000 m.為增加信號強度發(fā)射源長度設(shè)為100 m，發(fā)射波形為4階m序列偽隨機碼，碼元寬度為30 ms(圖18).我們選擇收發(fā)距分別為5000、5200、5500、5700 m和6000 m.

圖15 淺部高阻儲層對脈沖及階躍響應(yīng)影響 (a) 相關(guān)辨識脈沖響應(yīng)dex/dt；(b)積分階躍響應(yīng)ex.Fig.15 Effects of shallow hydrocarbon reservoir on impulse and step response (a) Impulse response by correlation identification; (b) Step response obtained by integration.

圖16 不同收發(fā)距反演結(jié)果Fig.16 Inversion results for different offsets

圖17 深部高阻儲集層模型Fig.17 Model for a deep hydrocarbon reservoir

圖18 4階m序列發(fā)射波形Fig.18 4-order m-sequence transmitting waveform

圖19是收發(fā)距為5000 m的m序列電場ex響應(yīng).其中黑色實線為不含有高阻層的均勻半空間響應(yīng)，灰色實線為含有高阻儲層的電磁響應(yīng).與圖14相比，由于高阻層埋深加大，收發(fā)距增加，高阻儲層的出現(xiàn)對電場ex的影響減弱，但在on-time和off-time與均勻半空間模型均有一定差異.在正演數(shù)據(jù)中加入高斯噪聲，信噪比為65 dB，然后利用與模型一相同的數(shù)據(jù)處理方法進行脈沖響應(yīng)辨識，并積分計算階躍響應(yīng).

圖20為深部高阻儲層存在與否對脈沖和階躍響應(yīng)影響對比圖.圖20a和圖20b分別是利用相關(guān)辨識技術(shù)獲得的脈沖響應(yīng)以及積分計算的階躍響應(yīng).其中黑色和灰色實線分別為不含高阻儲層的均勻半空間響應(yīng)以及含有高阻儲層的響應(yīng)曲線.與淺部高阻層的影響(圖15)相比，深部高阻層的出現(xiàn)同樣使得脈沖響應(yīng)產(chǎn)生較高的峰值，峰值出現(xiàn)時刻發(fā)生偏移，但幅度明顯減弱，且寬度變大.對階躍響應(yīng)影響的幅度也明顯減弱.

利用Occam方法對不同收發(fā)距的共中心點道集(CMP)數(shù)據(jù)進行反演，可以獲得多組反演結(jié)果(圖21).圖中黑色細實線、灰色細實線、黑色粗實線、黑色虛線和灰色粗實線分別是收發(fā)距為5000、5200、5500、5700 m和6000 m的反演結(jié)果，灰色虛線是對以上5個收發(fā)距數(shù)據(jù)同時反演結(jié)果.由圖可以看出，高阻層上下層的電阻率得到很好的反映，高阻層的深度和厚度也得到很好的反映.與淺部高阻層的反演結(jié)果相似(參見圖16)，利用5個收發(fā)距數(shù)據(jù)同時處理可以獲得更加理想的反演結(jié)果.因此，利用MTEM裝置同時測量多道數(shù)據(jù)探測深部高阻儲層同樣十分有效.

圖19 m序列電磁響應(yīng)Fig.19 EM responses for m-sequence

圖20 深部高阻儲層對脈沖及階躍響應(yīng)影響 (a)相關(guān)辨識脈沖響應(yīng)dex/dt；(b)積分階躍響應(yīng)ex.Fig.20 Effects of deep hydrocarbon reservoir on impulse and step responses (a) Impulse response by correlation identification; (b) Step response obtained by integration.

圖21 不同收發(fā)距反演結(jié)果Fig.21 Inversion results for different offsets

7 結(jié)論

多通道瞬變電磁方法采用m序列發(fā)射波形，同時觀測同線電場分量多道數(shù)據(jù)，在油氣儲層勘探中可取得較好的應(yīng)用效果.頻譜分析結(jié)果顯示m序列具有較寬的頻帶，譜線之間等間距，而且包含直流成分.通過改變碼元寬度和階數(shù)可以調(diào)整頻帶寬度以及能量分布，從而達到提高信噪比和分辨率的目的.根據(jù)實際勘探目標(biāo)合理調(diào)整m序列的相關(guān)參數(shù)，可以獲得理想的勘探效果.本文利用方波響應(yīng)移位疊加以及電流對時間導(dǎo)數(shù)與階躍響應(yīng)褶積的方法實現(xiàn)了對理論m序列和實際發(fā)射波形電磁響應(yīng)的全時模擬.采用相關(guān)辨識技術(shù)去除波形，計算脈沖響應(yīng).對其進行積分可以獲得階躍響應(yīng).m序列尖銳的相關(guān)特性，可以有效壓制噪聲干擾.直接對處理得到的階躍響應(yīng)數(shù)據(jù)進行反演，大大提高了工作效率.本文給出的詳細m序列發(fā)射波形正演模擬以及后續(xù)數(shù)據(jù)處理流程同樣適用于2D和3D模型，為MTEM電磁勘查系統(tǒng)研發(fā)提供理論保障，必將推動該項技術(shù)的發(fā)展和應(yīng)用.

致謝 作者向中國科學(xué)院電子研究所武欣及本電磁研究團隊的其他相關(guān)人員在文章準(zhǔn)備過程中提供的幫助致謝.特別向?qū)徃迦撕途庉媽Ρ疚奶岢龅慕ㄔO(shè)性意見表示感謝.

Chen S D, Lin J, Zhang S. 2012. Effect of transmitter current waveform on TEM response.ChineseJ.Geophys. (in Chinese), 55(2): 709-716, doi: 10.6038/j.issn.0001-5733.2012.02.035.Chen X D, Zhang Y, Zhang J, et al. 2012. The applications of HTc SQUID magnetometer to TEM.ChineseJ.Geophys. (in Chinese), 55(2): 702-708, doi: 10.6038/j. issn.0001-5733.2012.02.034.Chun S H, Chen R J, Geng M H. 2014. Review of the pseudo-random m sequence and its application in electrical prospecting of exploration geophysics.ProgressinGeophysics(in Chinese), 29(1): 439-446, doi: 10.6038/pg20140164.Constable S C, Parker R L, Constable C G. 1987. Occam′s inversion: A practical algorithm for generating smooth models from electromagnetic sounding data.Geophysics, 52(3): 289-300.

Guan S S, Lin J, Ji Y J, et al. 2012. Analysis of excitation signal airborne transient effect on grounded electrical-source electromagnetic response.ChineseJournalofRadioScience(in Chinese), 27(4): 766-772.He J S. 1998. 2nPseudo-random signal and its application. ∥ Annual of the Chinese Geophysical Society (in Chinese), Xi′an: Xi′an Map Press.

He J S, Tong T G, Liu J X. 2009. Mathematical analysis and realization ofansequence pseudo-random multi-frequencies signal.JournalofCentralSouthUniversity(ScienceandTechnology) (in Chinese), 40(6): 1666-1671.

Li M, Wei W B, Luo W B, et al. 2012. Time-domain spectral induced polarization based on pseudo-random sequence.PureandAppliedGeophysics, 170(12): 2257-2262, doi: 10.1007/ s00024-012-0624-z.

Li X, Xue G Q, Liu Y A, et al. 2012. A research on TEM imaging method based on synthetic-aperture technology.ChineseJ.Geophys. (in Chinese), 55(1): 333-340, doi: 10.6038/j.issn.0001-5733.2012.01.034.

Lin Z H, Chen S Y, Wang Y Y. 2009. m sequence and its application in communication.ModernElectronicsTechnique(in Chinese), 32(9): 49-51, doi: 10.3969/j.issn.1004-373X.2009.09.014.

Liu G M. 1998. Effect of transmitter current waveform on airborne TEM response.ExplorationGeophysics, 29(2): 35-41.

Pei Y F, Yin C C, Liu Y H, et al. 2014. Calculation and application of B-field for time-domain airborne EM.ProgressinGeophysics(in Chinese), 29(5): 2191-2196, doi: 10.6038/pg20140530.Smith R S, Lee T J. 2004. Asymptotic expansions for the calculation of the transient electromagnetic fields induced by a vertical magnetic dipole source above a conductive halfspace.PureandAppliedGeophysics, 161(2): 385-397.

Tang J T, Li F, Luo W B. 2007. Electrical fine-exploration transmitter design based on invert-repeated m-sequence.ProgressinGeophysics(in Chinese), 22(3): 994-1000.Tang J T, Luo W B. 2008. Pseudo-random electromagnetic exploration based on Invert-Repeated m-Sequence correlation identification.ChineseJ.Geophys. (in Chinese), 51(4): 1226-1233, doi: 10.3321/j.issn:0001-5733.2008.04.034.Wright D, Ziolkowski A, Hobbs B. 2002. Hydrocarbon detection and monitoring with a multicomponent transient electromagnetic (MTEM) survey.TheLeadingEdge, 21(9): 852-864, doi: 10.1190/1.2144389.Wright D, Ziolkowski A. 2007. Suppression of noise in MTEM data. ∥ SEG Technical Program Expanded Abstracts, 549-553.

Wiener N. 1949. Extrapolation, Interpolation, and Smoothing of Stationary Time Series. New York: Wiley.

Wu M J, Du T C. 2002. The pseudonoise sequence code and the birth for computer.JournalofLiaoningTechnicalUniversity(NaturalScience) (in Chinese), 21(2): 203-206.

Xie X K, Yang X R. 2010. Analysis of anti-jamming to Pseudo-Random code fuse signal.ElectronicInformationWarfareTechnology(in Chinese), 25(3): 48-53.

Xue G Q, Chen W Y, Zhou N N, et al. 2013. Short-offset TEM technique with a grounded wire source for deep sounding.ChineseJ.Geophys. (in Chinese), 56(1): 255-261, doi: 10.6038/cjg20130126.

Yin C C, Huang W, Ben F. 2013. The full-time electromagnetic modeling for time-domain airborne electromagnetic system.ChineseJ.Geophys. (in Chinese), 56(9): 3153-3162, doi: 10.6038/cig20130928.

Ziolkowski A, Hobbs B A, Wright D. 2007. Multitransient electromagnetic demonstration survey in France.Geophysics, 72(4): F197-F209, doi: 10.1190/1.2735802.

Ziolkowski A, Parr R, Wright D, et al. 2010. Multi-transient electromagnetic repeatability experiment over the North Sea Harding field.GeophysicalProspecting, 58(6): 1159-1176, doi: 10.1111/j.1365-2478.2010.00882.x.

Ziolkowski A, Wright D, Mattsson J. 2011. Comparison of PRBS and square-wave transient CSEM data over Peon gas discovery, Norway. ∥ SEG Technical Program Expanded Abstracts, 583-588, doi: 10.1190/1.3628149.

附中文參考文獻

陳曙東, 林君, 張爽. 2012. 發(fā)射電流波形對瞬變電磁響應(yīng)的影響. 地球物理學(xué)報, 55(2): 709-716, doi: 10.6038/j.issn.0001-5733.2012 02.035.

陳曉東, 趙毅, 張杰等. 2012. 高溫超導(dǎo)磁強計在瞬變電磁法中的應(yīng)用. 地球物理學(xué)報, 55(2): 702-708, doi: 10.6038/j.issn.0001-5733.2012.02.034.

淳少恒, 陳儒軍, 耿明會. 2014. 偽隨機m序列及其在電法勘探中的應(yīng)用進展. 地球物理學(xué)進展, 29(1): 439-446, doi: 10.6038/pg20140164.

關(guān)珊珊, 林君, 稽艷鞠等. 2012. 激勵信號對地-空瞬變電磁響應(yīng)的影響分析. 電波科學(xué)學(xué)報, 27(4): 766-772.

何繼善. 1998. 2n系列偽隨機信號及應(yīng)用. ∥ 中國地球物理學(xué)會年刊. 西安: 西安地圖出版社.

何繼善, 佟鐵鋼, 柳建新. 2009.an序列偽隨機多頻信號數(shù)學(xué)分析及實現(xiàn). 中南大學(xué)學(xué)報(自然科學(xué)版), 40(6): 1666-1671.

李貅, 薛國強, 劉銀愛等. 2012. 瞬變電磁合成孔徑成像方法研究. 地球物理學(xué)報, 55(1): 333-340, doi: 10.6038/j.issn.0001-5733.2012.01.034.

林智慧, 陳綏陽, 王元一. 2009. m序列及其在通信中的應(yīng)用. 現(xiàn)代電子技術(shù), 32(9): 49-51, doi: 10.3969/j.issn.1004-373X.2009.09.014.

裴易峰, 殷長春, 劉云鶴等. 2014. 時間域航空電磁磁場計算與應(yīng)用. 地球物理學(xué)進展, 29(5): 2191-2196, doi: 10.6038/pg20140530.

湯井田, 李飛, 羅維斌. 2007. 基于逆重復(fù)m序列的精細探測電法發(fā)送機設(shè)計. 地球物理學(xué)進展, 22(3): 994-1000.

湯井田, 羅維斌. 2008. 基于相關(guān)辨識的逆重復(fù)m序列偽隨機電磁法. 地球物理學(xué)報, 51(4): 1226-1233, doi: 10.3321/j.issn:0001-5733.2008.04.034.

吳明捷, 杜天蒼. 2002. 偽隨機碼及計算機的產(chǎn)生. 遼寧工程技術(shù)大學(xué)學(xué)報(自然科學(xué)版), 21(2): 203-206.

謝雪康, 楊曉蓉. 2010. 復(fù)雜電磁環(huán)境下偽隨機信碼引信信號抗干擾性分析. 電子信息對抗技術(shù), 25(3): 48-53.

薛國強, 陳衛(wèi)營, 周楠楠等. 2013. 接地源瞬變電磁短偏移深部探測技術(shù). 地球物理學(xué)報, 56(1): 255-261, doi: 10.6038/cjg20130126.

殷長春, 黃威, 賁放. 2013. 時間域航空電磁系統(tǒng)瞬變?nèi)珪r響應(yīng)正演模擬. 地球物理學(xué)報, 56(9): 3153-3162, doi: 10.6038/cig20130928.

(本文編輯 何燕)

Multi-transient EM full-time forward modeling and inversion of m-sequences

QI Yan-Fu1, YIN Chang-Chun1*, WANG Ruo2, CAI Jing1

1CollegeofGeo-explorationSciencesandTechnology,JilinUniversity,Changchun130021,China2InstituteofGeologyandGeophysics,ChineseAcademyofSciences,Beijing100029,China

The traditional transient electromagnetic method is mainly used in mineral exploration. It is difficult to meet the demand for oil and gas exploration. The multi-transient electromagnetic (MTEM) system can solve this problem due to its grounded wire as transmitter and in-line electric field component, resulting in the underground resistive targets can be more easily detected. However, the domestic research for this method is still at the stage of theoretical derivation, lacking forward modeling and data processing. We work on full-time forward modelling and inversion for m-sequence MTEM method to provide theoretical basis for the on-going development of MTEM system as well as its data processing and interpretation.

10.6038/cjg20150731.

國家自然科學(xué)基金項目(41274121，41174111)、國家重大科研裝備研制項目(ZDYZ2012-1-03，ZDYZ2012-1-05-04)聯(lián)合資助.

齊彥福，男，1989年生，博士，主要從事地球物理電磁正反演理論和方法技術(shù)研究. E-mail：jdqiyanfu@126.com

*通訊作者 殷長春，男，1965年生，教授，國家“千人計劃”特聘專家，主要從事電磁勘探理論，特別是航空和海洋電磁方面的研究. E-mail: yinchangchun@jlu.edu.cn

10.6038/cjg20150731

P631

2014-11-24，2015-03-19收修定稿

齊彥福, 殷長春，王若等. 2015. 多通道瞬變電磁m序列全時正演模擬與反演.地球物理學(xué)報，58(7):2566-2577,

Qi Y F, Yin C C, Wang R, et al. 2015. Multi-transient EM full-time forward modeling and inversion of m-sequences.ChineseJ.Geophys. (in Chinese),58(7):2566-2577,doi:10.6038/cjg20150731.