奇異值分解提取陣列聲波時(shí)差的改進(jìn)方法

李鵬舉, 吳昀朔, 任 莉

(1.東北石油大學(xué) 地球科學(xué)學(xué)院, 黑龍江 大慶 163318; 2.東北石油大學(xué) 非常規(guī)油氣成藏與開發(fā)省部共建國(guó)家重點(diǎn)實(shí)驗(yàn)室培育基地， 黑龍江 大慶 163318; 3.中國(guó)石油集團(tuán)測(cè)井有限公司 大慶分公司解釋評(píng)價(jià)中心， 黑龍江 大慶 163000)

0 引 言

多極子陣列聲波測(cè)井是20世紀(jì)90年代和21世紀(jì)初出現(xiàn)的一種極具前景的地球物理測(cè)井技術(shù)，是目前應(yīng)用較為廣泛的現(xiàn)代測(cè)井方法之一[1-2]。多極子陣列聲波全波列波形曲線中包含縱波、橫波、斯通利波等豐富的曲線信息，這些信息可用于裂縫識(shí)別及評(píng)價(jià)、巖層滲透性評(píng)價(jià)、油氣層識(shí)別評(píng)價(jià)和地層各項(xiàng)異性分析等，因而很有必要提高提取時(shí)差信息的精度[3]。由于聲波全波列波形曲線主要是由縱波、橫波、偽瑞利波和斯通利波等四種組分波及隨機(jī)噪聲組成，各組分波不僅在時(shí)間上有分離現(xiàn)象，在頻率上也是分離的，且縱波與橫波都有它們的固有頻率，即各組分波都有其自身的相似性[4]。早期聲波測(cè)井?dāng)?shù)據(jù)處理方法主要是閾值法，這種方法的缺點(diǎn)是在有噪聲脈沖存在的情況下首波檢測(cè)容易出錯(cuò)[5-6]。后期Kimball與Marzetat提出了時(shí)差—時(shí)間相關(guān)分析法—STC方法[7-8]，該方法對(duì)速度接近的一些波組分的分辨力有限,并且該方法對(duì)波列中較弱的成分靈敏度較低[9]。因此，筆者采用奇異值分解法提取聲波時(shí)差，利用改進(jìn)后的奇異值分解法進(jìn)行陣列聲波數(shù)據(jù)處理，以提高時(shí)差提取的精度。

1 陣列聲波時(shí)差的提取原理

奇異值分解作為一種在工程理論中獲得廣泛應(yīng)用的數(shù)學(xué)方法，其具有很明顯的物理意義。奇異值分解方法在處理占據(jù)大量?jī)?chǔ)存空間并且十分復(fù)雜的矩陣時(shí)具有極大優(yōu)勢(shì)，用幾個(gè)結(jié)構(gòu)簡(jiǎn)單占據(jù)儲(chǔ)存量小的子矩陣相乘來對(duì)原始矩陣進(jìn)行表示，且對(duì)原始矩陣的重要的特性進(jìn)行詳細(xì)描述。

一組由n條波形曲線組成的全波列聲波測(cè)井記錄，每條波形曲線的采樣點(diǎn)數(shù)為m，將這組數(shù)據(jù)轉(zhuǎn)化為二維數(shù)組由矩陣X表示，其中各元素為xij，即：

(1)

矩陣X的奇異值分解等式可表示為

(2)

式中：U——方陣XXT的m×r階特征向量(ui)矩陣；

V——方陣XTX的n×r階特征向量(vi)矩陣；

Λ——對(duì)角線矩陣，矩陣中各個(gè)元素為矩陣X的奇異值σi(XXT或XTX的非負(fù)平方根)，各個(gè)元素按遞減順序，Λ=diag(σ1,σ2,…,σi);

r——矩陣X的秩。

矩陣XXT(或XTX)為協(xié)方差矩陣，由于協(xié)方差矩陣所具有的正定性質(zhì)，使矩陣XXT(或XTX)的特征值為實(shí)數(shù)且為非負(fù)數(shù)，因此，矩陣X的奇異值一定是存在的[10]。

(3)

式中：XP——當(dāng)?shù)赑項(xiàng)至最后一項(xiàng)奇異值σP+1，σP+2，…，σN取零值時(shí)的重構(gòu)矩陣，即XP=U·ΛPVT，其中ΛP=diag(σ1,σ2,…,σP,0,…,0)；

‖X‖F(xiàn)——矩陣X的F范數(shù)[11-12]。

2 改進(jìn)的陣列聲波時(shí)差提取方法

(4)

2.1 子空間類算法

通過奇異值分解算法得到的奇異值矩陣主對(duì)角線元素σk是按遞減順序排列的，奇異值矩陣主對(duì)角線元素σk從首項(xiàng)開始對(duì)相鄰2項(xiàng)依次進(jìn)行差分得到梯度序列為

B(k)=σk-σk+1，(k=1,2,…,M-5),

(5)

式中,B(k)——差分梯度序列第k個(gè)元素。

算法的基本原理為尋找梯度序列每一項(xiàng)元素與后5項(xiàng)元素之和的比值序列A(k)的最大值，并將最大值的序列號(hào)k作為空間分界點(diǎn)

(6)

根據(jù)性質(zhì)0≤A(k)≤1，當(dāng)A(k)最大時(shí)，其k值代入式(7)便可對(duì)信噪比進(jìn)行優(yōu)化，

(7)

式中：σSi——奇異值中有用信號(hào)部分；

σN——奇異值中無用信號(hào)部分。

2.2 加權(quán)函數(shù)法

在得到用于時(shí)差提取的信噪比αSNR后，將αSNR加一個(gè)加權(quán)函數(shù)，從而突出有用αSNR，提高時(shí)差提取的準(zhǔn)確性。

加權(quán)函數(shù)為

(8)

加權(quán)后的信噪比為

βGSNR=ρ·αSNR，

(9)

3 實(shí)際測(cè)井?dāng)?shù)據(jù)的處理

為驗(yàn)證改進(jìn)后的奇異值分解提取陣列聲波時(shí)差方法的實(shí)際應(yīng)用效果，應(yīng)用子空間類奇異值分解算法對(duì)瑪湖油田瑪十八井區(qū)某井XMAC陣列聲波數(shù)據(jù)進(jìn)行實(shí)際處理。

3.1 子空間類奇異值分解法時(shí)差提取步驟

步驟一利用拉格朗日插值法，在每?jī)蓚€(gè)記錄點(diǎn)之間插入五個(gè)點(diǎn)，將時(shí)間采樣間隔由原來12 μs變?yōu)? μs，大大提高了處理精度。

步驟二在實(shí)際測(cè)量過程中，測(cè)井聲系的選取由所測(cè)井段實(shí)際地層情況來決定。時(shí)窗寬度Tw的一般選取原則為：選取某種待檢測(cè)模式波的理論長(zhǎng)度作為時(shí)窗寬度。Tmin和Tmax分別代表任意一條波形曲線所含有有用信號(hào)的最小時(shí)間與最大時(shí)間；ΔTb、ΔTe、δT三個(gè)變量分別代表相鄰接收探頭間波形曲線的掃描時(shí)差ΔT的起始值、終止值和增量。

步驟三將一組聲波全波列波形數(shù)據(jù)以T=Tmin為時(shí)窗中心，以ΔTb為掃描時(shí)差ΔT的起始值，以Tw為時(shí)窗寬度形成二維數(shù)據(jù)X。對(duì)得到的矩陣XXT或XTX進(jìn)行奇異值σi的計(jì)算，隨后計(jì)算信噪比αSNR，并對(duì)其進(jìn)行優(yōu)化，優(yōu)化后的信噪比βGSNR為相鄰波形曲線相似性的度量。

步驟四令ΔT=ΔTb+δT·i(i=1,2,…)，重復(fù)步驟二,直到ΔT=ΔTe為止,得到一條以Tmin為參數(shù)的ΔT-GSNR曲線。

步驟五令T=Tmin+(Tw/2),重復(fù)步驟二、步驟三，得到以Tmin+(Tw/2)為參數(shù)的ΔT-βGSNR曲線，如圖1所示。

步驟六令T=Tmin+(Tw/2)·j(j=2,3,…)，重復(fù)步驟二、三、四，直到T=Tmax為止，將得到的曲線合成繪制成βGSNR與ΔT、T的函數(shù)關(guān)系，如圖2所示，函數(shù)關(guān)系圖的峰值對(duì)應(yīng)的時(shí)差值作為時(shí)差結(jié)果進(jìn)行輸出。

圖1 時(shí)差提取步驟Fig. 1 Extract interval transit time

圖2 函數(shù)關(guān)系Fig. 2 Function relationship

圖2中,橫坐標(biāo)代表時(shí)差ΔT，縱坐標(biāo)代表聲波到時(shí)t。從圖2中可以看出，縱波時(shí)差的大致分布范圍為131～262 μs/m；橫波時(shí)差的分布范圍大致在295～427 μs/m之間。

3.2 數(shù)據(jù)分析

通過對(duì)陣列聲波測(cè)井信息提取方法的研究和總結(jié)，采用改進(jìn)后的奇異值分解方法借助CIFLog陣列聲波處理模塊對(duì)實(shí)際XMAC陣列聲波數(shù)據(jù)進(jìn)行處理。圖3是瑪湖油田瑪十八井區(qū)某井在2 710至2 820 m井段模式波時(shí)差提取結(jié)果。從圖3中可以看出，基于改進(jìn)后的奇異值分解法所提取出的聲波縱波時(shí)差DTC-New曲線與CIFLog軟件STC方法提取出的縱波時(shí)差DTC曲線結(jié)果基本一致，沒有出現(xiàn)大規(guī)模的曲線異常現(xiàn)象，且曲線間相關(guān)性較好。

該井?dāng)?shù)據(jù)各模式波形相關(guān)系數(shù)較高，根據(jù)DTC-New曲線可知，2 710～2 830 m深度段時(shí)差值穩(wěn)定在187～243 μs/m之間，且曲線較為平滑，在2 719～2 732 m處以及2 811～2 818 m處時(shí)差曲線出現(xiàn)3～62 μs/m的波動(dòng)，是由于地層巖性突變以及軟地層造成的橫波缺失所導(dǎo)致的結(jié)果，對(duì)縱波時(shí)差提取產(chǎn)生一定影響。在模式波時(shí)差曲線道橫波時(shí)差曲線與偶極源橫波曲線基本重合，總體上，根據(jù)實(shí)際測(cè)井情況可得出結(jié)論，基于改進(jìn)后的奇異值分解法提取的陣列聲波測(cè)井縱波時(shí)差結(jié)果較好，精度較高。

4 結(jié)束語

在陣列聲波測(cè)井時(shí)差提取的方法基礎(chǔ)上，對(duì)常規(guī)的奇異值分解法加以改進(jìn)，通過子空間類算法和加權(quán)函數(shù)法，大大提高了聲波時(shí)差的提取精度，應(yīng)用子空間類奇異值分解法對(duì)瑪湖油田瑪十八井區(qū)某井XMAC測(cè)井?dāng)?shù)據(jù)進(jìn)行解釋處理。通過與STC方法時(shí)差提取結(jié)果對(duì)比分析可知，改進(jìn)后的奇異值分解法提取出的時(shí)差結(jié)果更接近地層真實(shí)數(shù)據(jù)。該方法能夠提取出可靠的聲波時(shí)差，且計(jì)算效率較高，表明該方法在大多數(shù)地層中具有一定的實(shí)用價(jià)值，為進(jìn)一步的陣列聲波測(cè)井資料處理提供了新的思路與方法。