高 怡,李東航,曹景致,張翔宇,羅 旋

(1.西安石油大學(xué) 電子工程學(xué)院,陜西 西安 710065; 2.中國石油集團測井有限公司 測井技術(shù)研究院,陜西 西安 710077; 3.中國石油集團測井有限公司 吉林分公司,吉林 松原 138000)

引 言

在鹽水泥漿和高阻地層環(huán)境條件下,側(cè)向測井在地層電阻率測量和油氣藏評價中發(fā)揮著至關(guān)重要的作用[1]。相比雙側(cè)向儀器,陣列側(cè)向縱向分辨率高,測量信息更為豐富,5種探測深度分布合理,能夠清晰描述從井眼到原狀地層電阻率的變化情況,可為準確識別油氣層提供更可靠的資料[2]。

想要實現(xiàn)5種探測深度的“陣列式”測井,儀器需同時產(chǎn)生5種不同頻率的正弦測井信號。而這也就意味著,該儀器的測井響應(yīng)信號是一個多頻正弦波混合信號。如何準確、快速地提取出每個所需頻率處信號的幅值,成為了最終準確計算地層視電阻率的關(guān)鍵。在此背景下相敏檢波技術(shù)得到廣泛的應(yīng)用。

目前,井下儀器相敏檢波由模擬乘法器和濾波器構(gòu)成的模擬相敏檢波向數(shù)字相敏檢波(Digital Phase-Sensitive Detection,DPSD)轉(zhuǎn)變,既減少了模擬電路模塊,又提高了運算速度和線性度[3-4]?？焖俑道锶~變換(Fast Fourier Transform,FFT)也是一種應(yīng)用廣泛的數(shù)字檢波方法,在通信、信號頻譜分析等領(lǐng)域發(fā)揮著重要作用[5]。

本文以陣列側(cè)向測井儀為背景,分別介紹了目前幾種主流混頻信號數(shù)字檢波方法在儀器中的應(yīng)用,并針對傳統(tǒng)基-2 FFT算法中的倒位序算法進行了優(yōu)化。通過公式推導(dǎo)與測試,對比了幾種算法在運算速度與計算精度上的差異。最后以一段標準井的電阻率實測曲線分析了DPSD與改進的FFT這兩種方法在儀器上的實際應(yīng)用效果。

1 陣列側(cè)向測井原理

為了在復(fù)雜測井環(huán)境以及地層電阻率與泥漿電阻率反差極大的情況下得到地層的真實電阻率信息,陣列側(cè)向測井儀器應(yīng)運而生[6]。陣列側(cè)向儀器采用多個頻率點測量,不同的頻率點對應(yīng)不同的電流模式,且不同頻率點之間互不干擾,所以同時可以完成所有頻率點獨立聚焦測量。

陣列側(cè)向測井儀的電極系結(jié)構(gòu)如圖1所示。A0為主電極,兩側(cè)的電極A1(A1′)～A6(A6′)在各種探測模式中分別作為屏蔽電極或回路電極,M0(M0′)～M5(M5′)為成對的監(jiān)控電極。

圖1 陣列側(cè)向測井儀電極系Fig.1 Electrode system of array lateral logging tool

5種探測模式(外加RAL0模式——泥漿電阻率測量模式)的工作原理如下[7]:

(1)RAL0模式:A0發(fā)射主流,電流返回到A1(A1′)～A6(A6′)。由于該模式主流不聚焦,大部分電流將在井眼泥漿中流動,所以RAL0測井曲線主要反映泥漿電阻率值。

(2)RAL1模式:A0發(fā)射主流,A1(A1′)發(fā)射屏蔽電流,所有電流返回到A2(A2′)～A6(A6′),聚焦條件:M0(M0′)與M1(M1′)等電位。由于返回電極距A0不遠,主流進入地層較淺的地方就散開了,所以該模式的探測深度較淺。

(3)RAL5模式:A0發(fā)射主流,電極A1(A1′)～A5(A5′)發(fā)射屏蔽電流,而所有電流返回到A6(A6′)電極,聚焦條件為:M0(M0′)與M1(M1′)、M2(M2′)與M3(M3′)、M4(M4′)與M5(M5′)、A3(A3′)與A4(A4′)、A4(A4′)與A5(A5′)兩兩等電位,該模式探測深度最深。

(4)其他模式的工作原理可由RAL1模式與RAL5模式類比得出。

2 測井響應(yīng)信號特點與信號采集系統(tǒng)

陣列側(cè)向測井中,數(shù)字檢波需要檢測的測井響應(yīng)信號是由5種頻率正弦波混合信號或單一頻率正弦波信號所組成?；祛l信號成分的頻率分別為65、130、195、260以及325 Hz,各頻率間隔65 Hz,且各正弦信號間均滿足三角正交性。

具體來說,陣列側(cè)向儀器需要采集的測井響應(yīng)信號包含如下3類:

(1)監(jiān)督電極M0、M1b間的電位差信號ΔVM0M1b:該信號為5種頻率正弦波混合波形,數(shù)字檢波分離出這5種單頻信號的幅度值后,用于主電流聚焦計算及屏流或返回電極等電位監(jiān)控;

(2)監(jiān)督電極M0的電位信號VM0:同樣為5種頻率正弦波混合波形,數(shù)字檢波分離出5種單頻信號幅值用于視電阻率的計算;

(3)主電極A0的電流信號IA0:單頻率正弦波,需檢測出其幅值用于視電阻率計算。

在獲得以上幾類信號幅值后可通過

(1)

計算地層視電阻率。其中:i=1,2,3,4,5,分別代表5種頻率;K為儀器系數(shù),與儀器外殼材質(zhì)、各電極尺寸及直徑有關(guān)。

儀器采回的測井響應(yīng)信號首先通過信號調(diào)理電路進行預(yù)處理,然后由A/D芯片進行模數(shù)轉(zhuǎn)換,FPGA器件主要用于時序邏輯和組合邏輯控制,數(shù)字信號最終送入DSP芯片中進行相應(yīng)的數(shù)字檢波計算,并與儀器遙傳進行通信,地層實時的電阻率將在上位機測井軟件中呈現(xiàn)出來。該采集系統(tǒng)結(jié)構(gòu)框架如圖2所示。

圖2 儀器測井響應(yīng)信號采集系統(tǒng)Fig.2 Logging response signal acquisition system of logging instrument

3 陣列側(cè)向中的數(shù)字檢波

3.1 數(shù)字相敏檢波

數(shù)字相敏檢波與早期的模擬相敏檢波相比,克服了傳統(tǒng)模擬檢波器中由于三極管管壓降引起檢波精度不高的弊病,充分利用了測量信號實部和虛部的信息,使得儀器計算的地層電阻率更加精確,進一步提高了儀器精度[8]。相敏檢波原理如圖3所示[9]。

圖3 相敏檢波原理示意圖Fig.3 Schematic diagram of phase-sensitive detection principle

假設(shè)有一個正余弦三角函數(shù)集合{1,cosωt,…,cosmωt,sinωt,…,sinnωt,…}在區(qū)間(t0,t0+T)內(nèi)組成正交函數(shù)集,其中t0為初始時刻,基頻周期T=2π/ω,則以下關(guān)系成立:

(2)

數(shù)字相敏檢波實際上是用一個和待測量信號頻率嚴格相同的標準參考信號乘以該待測信號,乘積結(jié)果由兩部分組成——頻率為待測量信號兩倍的正弦信號,以及大小與待測量信號成比例的直流分量。

設(shè)Acosα為輸入的待測量正弦信號,2cosβ為同頻率參考信號,其中α=2πFt+θ,β=2πFt,角頻率ω=2πF。兩者相乘有:

2Acosαcosβ=A[cos (α+β)+cos (α-β)]=A[cos (4πFt+θ)+cosθ]。

(3)

為了消除上式中的第一項交流項,將式(3)在tk,k=0,1,…,N-1處離散化并求和:

式中:Fs為采樣頻率。

當待測信號與標準參考信號嚴格同相時,即θ=0時,式(4)的第一項在一個或多個周期內(nèi)離散求和(相當于對三角函數(shù)求積分)出現(xiàn)正負相消,結(jié)果為0,第二項中的cosθ值為1,故整個式子的結(jié)果即為該信號的幅值A(chǔ),即“檢”出信號的幅值。

陣列側(cè)向測井儀中的待測信號可設(shè)為:

x(t)=A1cos 2πF1t+A2cos 2πF2t+A3cos 2πF3t+A4cos 2πF4t+A5cos 2πF5t。

(5)

式中:F1～F5為陣列側(cè)向儀器中的5種頻率信號,各信號對應(yīng)幅度為A1～A5。

對式(4)兩端乘以頻率為Fi的標準正弦信號cos 2πFit,并在時間Δt內(nèi)離散求和得到:

(6)

由式(2)正交函數(shù)性質(zhì),式(6)的結(jié)果最終為:

(7)

該結(jié)果即為某頻率處的信號幅值。

3.2 快速傅里葉變換及其優(yōu)化

對于N點序列x(n),其離散傅里葉變換(DFT)定義為[10]:

(8)

顯然,求取N點,X(k)需要N2次復(fù)數(shù)乘法及N(N-1)次復(fù)數(shù)加法。當N較大時,計算量非常龐大?？焖俑道锶~變換是DFT的快速算法,能極大減小復(fù)數(shù)乘加的計算量。

令N=2M,M為正整數(shù),令n=2r及n=2r+1,r=0,1,…,N/2-1,則式(8)可以按奇、偶數(shù)分成兩組,并且使用“折半定理”有:

(9)

以N/2作為中心點,設(shè)a為小于N/2的一個偶數(shù),A為a的逆序數(shù)(A≠a),且利用逆序數(shù)的對偶性質(zhì)得到下列關(guān)系:

A=rev(a),a=rev(A)。

(10)

其中,函數(shù)rev(·)為倒位序函數(shù)。

對于奇數(shù)則有如下關(guān)系:

(11)

更進一步,還需要如下關(guān)系:

(12)

利用如上的式(10)～式(12)即可直接計算得到輸出序列X(k)的倒置碼位。下面以N=16的情形為例來說明計算過程:中心點N/2=8,第一步由式(10)得到小于8的偶數(shù)所對應(yīng)的逆序數(shù);第二步,對于小于8的奇數(shù)只需在其前一個逆序數(shù)上加8,即式(11)的效果;第三步,根據(jù)對偶性,得到大于8的偶數(shù)逆序數(shù);最后利用式(12)計算得到其余逆序數(shù),計算過程如圖4所示。使用該優(yōu)化方法計算倒位序的計算量將遠小于雷德算法。

圖4 改進算法計算過程(以N=16為例)Fig.4 Calculation process diagram of improved algorithm

此外,FFT運算的結(jié)果還可得到各頻率測井響應(yīng)信號與標準正弦信號的相位差,這對于儀器模擬部分電路的設(shè)計具有一定的指導(dǎo)意義。

4 性能對比與測試

幅值提取的運行速度與計算精度是數(shù)字檢波方法的兩大性能指標。

使用C語言分別對上述DPSD、基-2 FFT以及改進的算法進行編程,針對N=221情形在同一計算機上進行運算速度測試。測試結(jié)果見表1。

表1 數(shù)字檢波算法運算速度測試結(jié)果Tab.1 Operation speed of digital detection algorithms

表1中的數(shù)據(jù)表明:DPSD算法運算速度與FFT相比具有明顯優(yōu)勢,運算速度為傳統(tǒng)基-2 FFT的7倍左右;改進倒位序后的FFT算法的運算速度較原來提高了近1倍,優(yōu)化效果明顯。值得指出的是,使用DPSD對多頻率混頻信號進行幅值提取時,需要采用多通道模式進行并行運算,這樣會使數(shù)據(jù)處理速度變緩,且會占用較多的芯片資源。

另外,為對比各檢波方法的幅值提取計算精度,基于TMS320F28335型DSP芯片進行編程、下載,搭配同一A/D芯片,分別使用DPSD與基-2 FFT對相同單一頻率信號進行幅值提取,改進算法主要影響FFT的計算速度,所以未在此處進行對照。分別設(shè)置待測信號幅值為500 mV、100 mV與10 mV,在每個幅值下重復(fù)進行5組實驗,DPSD與FFT的測試結(jié)果對比見表2,最后一列計算了5組數(shù)據(jù)的平均值與真值間的相對誤差,用于直觀對比計算精度。實驗中為控制變量,DSP主頻、采樣點數(shù)均保持一致。

表2 基于DSP的FFT和DPSD精度測試結(jié)果Tab.2 Accuracy test results of FFT and DPSD based on DSP

由表2數(shù)據(jù)對比可以得出結(jié)論:當所提供的待測信號與標準信號沒有相位差時,DPSD方法的幅值提取運算精度高于同點數(shù)的FFT算法,特別是對極小信號的幅值提取優(yōu)勢很大,所以DPSD方法適用于對小信號要求比較高的測井儀器。產(chǎn)生此現(xiàn)象的原因是DPSD原理上是同頻檢測,具有選頻能力,從而提高了抗干擾能力。但是當待測信號的相位與標準信號相位存在偏差時,DPSD將產(chǎn)生巨大誤差。圖5所示是DPSD與FFT對于500 mV信號的幅值提取計算精度隨信號相移大小的變化情況。

圖5 運算誤差與相移關(guān)系Fig.5 Relationship between error and phase shift

圖5中,隨著信號相移的不斷增大,DPSD的計算誤差將大幅上升,而FFT方法的誤差始終保持在1%以內(nèi)。其原因為:FFT有子載波和周圍頻點的幅值信息可以后期彌補和糾正,而DPSD只針對某單頻來計算,會把其余頻率信息丟失掉。

由于陣列側(cè)向儀器的模擬電路部分會導(dǎo)致測井響應(yīng)信號產(chǎn)生相移,在DPSD與FFT的運算速度都能滿足儀器數(shù)據(jù)上傳速率的情況下,FFT能保證在發(fā)生相移時測井信號的幅值提取運算精度,在陣列側(cè)向測井儀中得到了廣泛應(yīng)用。

5 陣列側(cè)向儀器應(yīng)用實例

為對比DPSD與改進的FFT在陣列側(cè)向儀器上的實際應(yīng)用效果,分別將2種算法下載至儀器主控板DSP芯片中并裝入儀器,在X標準井進行下井實測。

該標準井井深為612 m,井內(nèi)有直徑為8 inch(203.2 mm)鋼套管,內(nèi)徑200 mm,下井測試的陣列側(cè)向儀器外直徑55 mm,截取的儀器所測該井383～405 m段的電阻率曲線如圖6所示,該段地質(zhì)成分為粉質(zhì)黏土及中細沙互層。由于儀器直徑較小,遠小于井眼內(nèi)徑,所以RAL1電阻率曲線受泥漿影響較大,并且 RAL0模式僅僅反映了泥漿電阻率值,參考價值不高,故圖中刪除了RAL0、RAL1,只給出了其余4條曲線。

圖6 X標準井電阻率實測曲線Fig.6 Measured resistivity curve of X standard well

圖6中左右兩道曲線是分別使用2種數(shù)字檢波對同一井段的電阻率實測曲線,左側(cè)為FFT,右側(cè)為DPSD。圖中的紅色點劃線是使用測井軟件選取的RAL5曲線最高點,旁邊顯示的是該點所處的深度與電阻率值。

由該段測井曲線可以看出:與DPSD相比,應(yīng)用FFT方法所測得的4條電阻率曲線分布更加合理,測井曲線的一致性更好,測量的數(shù)值符合該井段地質(zhì)結(jié)構(gòu);FFT方法在393 m井深處能測出視電阻率值高達45 000 Ω·m的高阻地層,儀器測量精度更高,而DPSD在此井深處的測量值僅為18 000 Ω·m左右,測量偏差較大。另外,在393 m處曲線出現(xiàn)負差異現(xiàn)象,可能與絕緣棒影響、圍巖影響、地層的各向異性等因素有關(guān)。

6 結(jié) 論

(1)DPSD算法運算速度與FFT相比具有明顯優(yōu)勢,運算速度大約為傳統(tǒng)基-2 FFT的7倍左右,改進倒位序后的基-2 FFT算法的運算速度較原來提高了近1倍;

(2)當待測信號與標準信號無相位差時,DPSD算法的幅值提取計算精度略高于FFT,隨著相移增大,DPSD的計算誤差會大幅上升,而FFT的誤差則不受相位偏移的影響;

(3)在X標準井測井實驗中,儀器在FFT算法下所測的電阻率曲線一致性更好,測量精度更高,能測出視電阻率值高達45 000 Ω·m的高阻地層,驗證了FFT更適合陣列側(cè)向儀器的應(yīng)用,這將為準確識別油氣層提供更加可靠的資料。