李擁軍,宋 煒,唐傳章,史應(yīng)龍,王澤丹,陳樹(shù)光,劉 靜,王 標(biāo)
(1.中國(guó)石油天然氣股份有限公司華北油田巴彥勘探開(kāi)發(fā)分公司,內(nèi)蒙古巴彥淖爾750306;2.中國(guó)石油大學(xué)(北京)地球物理學(xué)院,北京102249;3.中國(guó)石油天然氣股份有限公司華北油田分公司勘探事業(yè)部,河北任丘062550;4.中國(guó)石油天然氣股份有限公司華北油田分公司勘探開(kāi)發(fā)研究院,河北任丘062550)
FUTTERMAN[1]在介質(zhì)對(duì)地震波吸收衰減效應(yīng)的研究中指出,巖石的吸收衰減特征是地層的基本屬性之一,因此地層的吸收衰減參數(shù)的精確求取問(wèn)題一直是研究的熱點(diǎn)。而與深部地層相比,近地表地層對(duì)地震波的吸收衰減作用更為強(qiáng)烈,但在處理過(guò)程中,很難準(zhǔn)確地求取地層吸收衰減參數(shù),尤其是近地表Q值參數(shù)。目前的研究已經(jīng)給出了多種Q值估計(jì)方法,且各種方法都有優(yōu)缺點(diǎn)和相應(yīng)的應(yīng)用范圍?;跀?shù)據(jù)來(lái)源的不同,可將已有的方法分為VSP數(shù)據(jù)和地面反射地震數(shù)據(jù)兩大類,而微測(cè)井可以歸為VSP數(shù)據(jù)類?;赩SP數(shù)據(jù)類的算法主要包括譜比法、質(zhì)心頻移法及峰值頻移法等。BATH[2]提出了譜比法,即提取兩個(gè)時(shí)間或深度處的子波頻譜,使用兩振幅譜譜比的斜率來(lái)計(jì)算Q值。GLADWIN等[3]基于衰減時(shí)脈沖寬度增加的現(xiàn)象提出了上升時(shí)間原理,KJARTANSSON[4]依據(jù)上升時(shí)間原理提出了采用上升時(shí)間法求Q值的方法。JANNSEN等[5]提出子波模擬法,通過(guò)改變Q值進(jìn)行正演,選擇正演結(jié)果與實(shí)際信號(hào)最接近的模擬Q值作為近似,之后依照近似原理又出現(xiàn)了頻率域近似的頻譜模擬法。TONN[6]使用合成VSP資料對(duì)不同Q值估計(jì)方法(包括譜比法、解析信號(hào)法、振幅衰減法、子波模擬法、擬合技術(shù)等)進(jìn)行了對(duì)比分析,認(rèn)為大多數(shù)方法求取的Q值精度在很大程度上取決于地震資料的質(zhì)量,且方法大都建立在特定假設(shè)的條件下,實(shí)際應(yīng)用時(shí)缺乏普適性。
近年來(lái),諸多學(xué)者在前人研究的基礎(chǔ)上提出了一些改進(jìn)算法,其中,趙秋芳等[7]梳理了現(xiàn)有近地表品質(zhì)因子估算方法,將其初步劃分為巖石樣本測(cè)試Q值估算和地層原位測(cè)量Q值估算兩大類;金子奇等[8]采用對(duì)數(shù)譜比反演方法,利用多道反射波信息同時(shí)反演多道衰減旅行時(shí),從而避免了傳統(tǒng)方法中人為選擇最優(yōu)頻帶帶來(lái)的誤差;陳文爽等[9]、王小杰等[10]和王德利等[11]基于S變換的時(shí)頻譜作為譜比法的輸入消除傳統(tǒng)譜比法在計(jì)算地層吸收參數(shù)時(shí)的平均效應(yīng);張繁昌等[12]采用自適應(yīng)分解的子波通過(guò)譜比法來(lái)分析薄層干涉問(wèn)題;李偉娜等[13]借鑒微測(cè)井分層速度回歸分析思想,提出了一種雙線性回歸穩(wěn)定Q值估計(jì)算法提高了精度;王宗俊[14]利用譜模擬獲得的子波譜約束質(zhì)心頻移法進(jìn)行Q值估算;張立彬等[15]引入積分中值參變量避免了對(duì)質(zhì)心頻移法中震源子波的先驗(yàn)假設(shè);高靜懷等[16]利用特征結(jié)構(gòu)法改善地震子波峰值頻率提取精度提高Q值估算精度;馮瑋等[17]采用時(shí)間域子波匹配法通過(guò)子波褶積矩陣引入衰減響應(yīng),從根本上避免譜估計(jì)的誤差;郭銳等[18]引入Capon2DQ值估計(jì)新方法,對(duì)其求解過(guò)程進(jìn)行了加權(quán)改進(jìn);宋吉杰等[19]提出一種基于信息融合的近地表介質(zhì)Q值估計(jì)方法以及穩(wěn)定的反Q濾波處理技術(shù),提高了地震資料的分辨率;王靜等[20]利用零井源距VSP資料,采用改進(jìn)的頻譜擬合法進(jìn)行穩(wěn)定的深層Q值估算;蘇勤等[21]基于地表一致性原理,引入表層相對(duì)衰減系數(shù)的概念,在共炮檢域迭代計(jì)算相對(duì)衰減系數(shù)并實(shí)現(xiàn)穩(wěn)定的表層Q值空變補(bǔ)償。上述方法重點(diǎn)是在Q值求取算法中進(jìn)行了各種改進(jìn),本文重點(diǎn)則是要在微測(cè)井地震數(shù)據(jù)初至波分離、走時(shí)自動(dòng)提取和頻譜求取方面進(jìn)行改進(jìn)。MALLAT[22]首次提出了一種使用Gabor小波作為時(shí)頻原子庫(kù)的匹配追蹤方法。LIU等[23]進(jìn)一步提出了使用Ricker子波庫(kù)和Morlet子波庫(kù)的快速匹配追蹤分解方法,并在Morlet原子庫(kù)中將子波尺度參數(shù)作為常數(shù)。WANG[24]將尺度參數(shù)改為可變參數(shù),并對(duì)Morlet原子庫(kù)進(jìn)行了改進(jìn),使殘差下降速度更快,完備庫(kù)中的原子數(shù)更多。SONG[25]提出了一種基于快速匹配追蹤的近地表Q值估計(jì)方法,提高了Q值估計(jì)的精度和抗噪能力,但使用零相位Ricker小波原子在匹配追蹤提取初至?xí)r會(huì)產(chǎn)生一定誤差。在無(wú)噪情況下,常規(guī)的譜比法是可靠的方法,當(dāng)信噪比較高時(shí),譜比法依賴于諸多因素,如時(shí)窗長(zhǎng)度、時(shí)窗形狀、求取斜率的起止頻率等,當(dāng)信噪比較低時(shí),譜比法可靠性差,因此,需要深入研究譜比法解的穩(wěn)定性。本文基于非零相位雷克子波庫(kù)的復(fù)數(shù)域快速匹配追蹤算法,實(shí)現(xiàn)微測(cè)井資料初至波分離、走時(shí)自動(dòng)拾取和頻譜求取,并通過(guò)引入整形正則化算子和最小平方反演方法求取譜比值,提出了一種新的近地表Q值估計(jì)方法。
微測(cè)井觀測(cè)點(diǎn)在整個(gè)地震工區(qū)比較稀疏,使得利用稀疏采樣點(diǎn)信息構(gòu)建三維近地表模型較為困難,傳統(tǒng)的雙線性、雙三次樣條及克里金插值建模的方法很難獲得較好的建模效果。近年來(lái),機(jī)器學(xué)習(xí)和深度學(xué)習(xí)被廣泛應(yīng)用于多個(gè)領(lǐng)域,在石油勘探、開(kāi)發(fā)領(lǐng)域也取得較好的應(yīng)用效果。R?TH等[27]通過(guò)神經(jīng)網(wǎng)絡(luò)反演了地下速度,證明了機(jī)器學(xué)習(xí)方法能夠解決反演問(wèn)題。MOHAMED等[28]使用神經(jīng)網(wǎng)絡(luò)完成了地震數(shù)據(jù)的疊后反演。KORJANI等[29]通過(guò)深度學(xué)習(xí)算法預(yù)測(cè)了地下巖石的物理參數(shù)。DAHLKE等[30]提出使用深度學(xué)習(xí)方法直接根據(jù)原始地震記錄進(jìn)行斷層預(yù)測(cè),在模型數(shù)據(jù)中實(shí)驗(yàn)準(zhǔn)確性較高?;谏疃葘W(xué)習(xí)的多元非線性回歸方法,利用工區(qū)范圍內(nèi)和周邊工區(qū)大量的微測(cè)井?dāng)?shù)據(jù)作為樣本進(jìn)行深度學(xué)習(xí),在一個(gè)近地表?xiàng)l件相對(duì)穩(wěn)定的地區(qū),可以利用大量的已知樣本數(shù)據(jù)建立統(tǒng)一的多元非線性回歸算子作為初始回歸模型,再通過(guò)遷移學(xué)習(xí)方法,將目標(biāo)工區(qū)內(nèi)新獲得的微測(cè)井?dāng)?shù)據(jù)作為新的學(xué)習(xí)標(biāo)簽對(duì)多元非線性回歸進(jìn)行調(diào)優(yōu),不斷優(yōu)化多元非線性回歸算子,以保證近地表模型的建模精度。本文將工區(qū)內(nèi)不同位置的微測(cè)井?dāng)?shù)據(jù)估計(jì)得到的Q值函數(shù)作為已知樣本標(biāo)簽,并通過(guò)訓(xùn)練深度神經(jīng)網(wǎng)絡(luò)來(lái)構(gòu)建多元非線性回歸算子,實(shí)現(xiàn)了三維近地表品質(zhì)因子建模,并將建立的三維近地表模型用于地震資料處理,以驗(yàn)證方法的有效性。
常規(guī)的匹配追蹤算法要求在每次匹配中掃描全部的原子庫(kù),但原子庫(kù)是過(guò)完備的,掃描過(guò)程計(jì)算量較大,計(jì)算速度較慢。為了減少匹配中的掃描工作量,可以使用快速?gòu)?fù)數(shù)域匹配追蹤(Complex domain fast matching pursuit decomposition,CFMP)方法來(lái)提高計(jì)算效率。CFMP方法需要根據(jù)實(shí)地震記錄獲取復(fù)地震記錄,從復(fù)信號(hào)中獲取先驗(yàn)信息(如振幅包絡(luò)最大值處的時(shí)間、瞬時(shí)頻率以及瞬時(shí)相位),然后在先驗(yàn)信息的約束下掃描一小部分的原子庫(kù)得到最佳匹配,完成復(fù)信號(hào)的重構(gòu)后再將結(jié)果返回實(shí)數(shù)域,從而在一定程度上減小匹配中的計(jì)算量,提高算法效率。匹配追蹤算法的基本思想是將信號(hào)投影到一系列時(shí)頻原子上,選取的時(shí)頻原子能很好地匹配信號(hào)的局部特征,利用這些時(shí)頻原子精確地表示原始信號(hào)。
F(t)=f(t)+i·H(f(t))
(1)
其中,f(t)為實(shí)信號(hào),虛部H(f(t))為實(shí)信號(hào)的Hilbert變換,其實(shí)質(zhì)是實(shí)部的90°相移。
同理,假設(shè)g(t)為時(shí)頻原子,H(g(t))表示時(shí)頻原子的Hilbert變換,則復(fù)時(shí)頻原子的表達(dá)式為:
G(t)=g(t)+i·H(g(t))
(2)
匹配追蹤從初始?xì)埩縍0F=F開(kāi)始,首先求取殘量的瞬時(shí)振幅包絡(luò)、瞬時(shí)頻率和瞬時(shí)相位,將振幅包絡(luò)最大值處的時(shí)間、瞬時(shí)頻率以及瞬時(shí)相位分別表示為tj,fj和φj。通過(guò)n次匹配計(jì)算,得到:
RnF=〈RnF,Gγn〉Gγn+Rn+1F
(3)
RnF為復(fù)殘量。其中,Gγn為復(fù)時(shí)頻原子,且‖Gγn‖=1,〈RnF,Gγn〉表示信號(hào)與殘差的內(nèi)積。很容易得出Gγn正交于Rn+1F,因此有:
‖RnF‖2=|〈RnF,Gγn〉|2+‖Rn+1F‖2
(4)
為了使‖Rn+1F‖最小,則所選擇的復(fù)子波Gγn∈G使得|〈RnF,Gγn〉|最大,即最佳原子的選擇條件為:
(5)
通過(guò)(5)式選擇的時(shí)頻原子能使每次迭代計(jì)算所得的殘量最小。當(dāng)信號(hào)的局部與匹配原子性質(zhì)完全相同或者完全匹配時(shí),復(fù)信號(hào)在原子上的投影值〈RnF,Gγn〉為一實(shí)數(shù),即有real(〈RnF,Gγn〉)=〈RnF,Gγn〉,故匹配條件改進(jìn)為:
(6)
通過(guò)(6)式可計(jì)算出最佳匹配子波原子,即計(jì)算出最佳的時(shí)移、主頻、相位等參數(shù)。為了使匹配效果最佳,本文在構(gòu)建時(shí)頻原子庫(kù)時(shí),對(duì)前面所求參數(shù)(時(shí)移、主頻、相位)的一定鄰域內(nèi)構(gòu)建時(shí)頻原子字典庫(kù)W,W={Wγ|‖Wγ‖=1},其中:指標(biāo)集γ是集合Γ中的元素,復(fù)時(shí)頻原子都經(jīng)過(guò)單位化。經(jīng)過(guò)匹配條件(6)式,得到首次匹配的結(jié)果為:
R0F=real(〈R0F,Wγ1〉)Wγ1+R1F
(7)
其中,R1F是首次迭代的復(fù)殘量,可直接用于下一次迭代計(jì)算,省去了Hilbert變換所帶來(lái)的計(jì)算量,計(jì)算過(guò)程與前一次相同。通過(guò)設(shè)置迭代次數(shù)或誤差閾值,復(fù)匹配追蹤最后的結(jié)果表示為:
(8)
由于信號(hào)中存在各種噪聲以及其它干擾因素,使得使用瞬時(shí)信息求得的匹配參數(shù)存在誤差,本文使用復(fù)數(shù)域的最小二乘法對(duì)求得的匹配參數(shù)進(jìn)行校正,使匹配參數(shù)所確定的復(fù)雷克子波能更接近地震信號(hào)的局部信息。設(shè)RnF=AjWγ+Rn+1F,其中Wγ為利用上面瞬時(shí)參數(shù)確定所得的復(fù)時(shí)頻原子,Aj=|Aj|·exp(iθj)表示對(duì)子波原子的校正項(xiàng),其中包括振幅校正|Aj|與相位校正θj。為使校正后的原子子波能最佳匹配地震信號(hào),故‖Rn+1F‖取最小值,使用復(fù)數(shù)域最小二乘法,求得:
(9)
通過(guò)(9)式求得的Aj得到振幅和相位校正項(xiàng)后的匹配子波原子變?yōu)镚γ=|Aj|Wγ(t),其中匹配參數(shù)為γ={uj,fj,φj},φj=αj+θj為校正之后的相位。由(9)式可計(jì)算出最佳匹配子波原子,即計(jì)算出最佳的時(shí)移、主頻、相位等參數(shù)。對(duì)復(fù)信號(hào)求取其先驗(yàn)信息,包括:振幅包絡(luò)、瞬時(shí)頻率和瞬時(shí)相位。復(fù)信號(hào)的振幅包絡(luò)、瞬時(shí)頻率和瞬時(shí)相位的計(jì)算公式分別為:
(10)
(11)
(12)
其中,f*(t)為信號(hào)的虛部。n次迭代后的振幅系數(shù)為:
(13)
另外,匹配追蹤采用非零相位雷克子波構(gòu)成時(shí)頻原子庫(kù),具體構(gòu)建非零相位雷克子波庫(kù)分兩步。首先給出零相位雷克子波W的解析表達(dá)式如下:
W=[1-2(πft)2]e-(πft)2
(14)
然后,對(duì)(14)式零相位雷克子波W進(jìn)行相位旋轉(zhuǎn),得到相位旋轉(zhuǎn)φ的雷克子波Wφ表達(dá)式為:
Wφ=Wcosφ+Imag(H(W))sinφ
(15)
其中,符號(hào)Imag表示取虛部。進(jìn)一步地,以振幅包絡(luò)最大值處的瞬時(shí)頻率和瞬時(shí)相位為基準(zhǔn),按一定的鄰域范圍將瞬時(shí)頻率和瞬時(shí)相位左右擴(kuò)展,再基于(14)式和(15)式,形成局部非零相位雷克子波時(shí)頻原子庫(kù)。
常規(guī)的對(duì)數(shù)譜比法就是兩式作比值并求對(duì)數(shù),即
(16)
式中:Δt為地震波在某地層中傳播的雙程旅行時(shí);A1(f)為衰減后的子波頻譜值;A2(f)為衰減前的子波頻譜值;Q為地層的品質(zhì)因子;c為常量。顯然,譜比對(duì)數(shù)值為頻率f的線性函數(shù),其斜率為k=π(Δt/Q),通過(guò)線性回歸求出斜率k值即可確定Q值。其中關(guān)鍵是求取穩(wěn)定的對(duì)數(shù)譜比值。
為了提高譜比計(jì)算的穩(wěn)定性和抗噪能力,可采用優(yōu)化穩(wěn)定反演算法來(lái)求取譜比參數(shù)。匹配追蹤方法可用于信號(hào)的分解,能夠提取不同時(shí)刻的子波,并獲得相應(yīng)的幅度譜。子波的譜比表達(dá)式可寫為:
(17)
式中:AiW(ti,f)是在ti處匹配得到的最佳原子振幅譜。當(dāng)分母數(shù)值很小時(shí),除法運(yùn)算的結(jié)果會(huì)出現(xiàn)不穩(wěn)定,簡(jiǎn)單的方法是在分母上加一個(gè)白噪系數(shù)(ξ)來(lái)防止不穩(wěn)定,即:
(18)
其中,ξ的選取和信號(hào)的信噪比有關(guān)。FOMEL[26]利用整形正則化來(lái)穩(wěn)定反演問(wèn)題,通過(guò)引入整形正則化的方法來(lái)獲得較為穩(wěn)定的譜比結(jié)果。通過(guò)引入正則化算子,將(18)式轉(zhuǎn)化成最小平方問(wèn)題:
(19)
式中:D是Tikhonov正則化算子;ti,ti+1分別為微測(cè)井不同深度位置的初至波旅行時(shí)。采用共軛梯度迭代來(lái)完成正則化反演,譜比的線性模型可以寫成:
m=L-1d
(20)
式中:m是O(f)的矩陣形式;d是Ai+1W(ti+1,f)的矩陣形式;L是AiW(ti,f)的矩陣形式。模型的理論解寫為:
(21)
引入整形正則化算子S=(I+ξ2DTD)-1,則正則化的理論解可以寫為:
(22)
其中,λ為阻尼調(diào)整因子。
神經(jīng)網(wǎng)絡(luò)是學(xué)習(xí)算法的一種,由多層互連節(jié)點(diǎn)組成“網(wǎng)絡(luò)”。受動(dòng)物神經(jīng)系統(tǒng)的啟發(fā),神經(jīng)網(wǎng)絡(luò)的網(wǎng)絡(luò)節(jié)點(diǎn)類似于神經(jīng)元,邊緣類似于突觸,如圖1所示。
圖1 生物神經(jīng)元(a)與人工神經(jīng)網(wǎng)絡(luò)(b)示意
在每個(gè)邊緣都有一個(gè)相關(guān)權(quán)重,網(wǎng)絡(luò)定義了輸入數(shù)據(jù)從網(wǎng)絡(luò)的輸入層傳遞到輸出層的計(jì)算規(guī)則。用神經(jīng)網(wǎng)絡(luò)的網(wǎng)絡(luò)函數(shù)表征輸入和輸出層之間的關(guān)系,由權(quán)重進(jìn)行參數(shù)化。適當(dāng)定義網(wǎng)絡(luò)函數(shù)后,可以通過(guò)最小化網(wǎng)絡(luò)權(quán)重的代價(jià)函數(shù)來(lái)執(zhí)行學(xué)習(xí)任務(wù)。多層的神經(jīng)網(wǎng)絡(luò)可以進(jìn)行特征學(xué)習(xí),網(wǎng)絡(luò)在隱藏層學(xué)習(xí)輸入,隨后在輸出層進(jìn)行分類或回歸。
微測(cè)井資料是精確求取近地表模型的關(guān)鍵,但是通常情況下由于成本的原因,微測(cè)井的空間采樣密度一般是在1~4km2內(nèi)布設(shè)一口井,在一塊100km2的三維工區(qū)內(nèi),通常微測(cè)井點(diǎn)有50~100個(gè),樣本點(diǎn)過(guò)少很難得到好的深度學(xué)習(xí)模型,因此在開(kāi)展多元非線性三維回歸模型訓(xùn)練過(guò)程中,盡可能多地收集工區(qū)附近的微測(cè)井?dāng)?shù)據(jù),作為初始模型訓(xùn)練的樣本,在此基礎(chǔ)上,基于遷移學(xué)習(xí)原理,利用本工區(qū)的樣本點(diǎn)深化調(diào)優(yōu),提高模型預(yù)測(cè)的精度。深度學(xué)習(xí)訓(xùn)練分為4個(gè)步驟:數(shù)據(jù)集整理、定義模型、訓(xùn)練/學(xué)習(xí)和預(yù)測(cè)/評(píng)估。數(shù)據(jù)集是本工區(qū)用前述譜比法估算的每個(gè)觀測(cè)點(diǎn)的Q值函數(shù),由測(cè)點(diǎn)坐標(biāo)(X,Y),測(cè)點(diǎn)深度(Depth)和Q值構(gòu)成,其中X、Y和深度構(gòu)成訓(xùn)練集的輸入X_train,對(duì)應(yīng)的Q值構(gòu)成訓(xùn)練集的標(biāo)簽label_train。我們需要對(duì)輸入數(shù)據(jù)進(jìn)行標(biāo)準(zhǔn)化,使其滿足正態(tài)分布,開(kāi)源的機(jī)器學(xué)習(xí)算法包Scikit-learn中包含的StandardScalar函數(shù)可以對(duì)任何輸入數(shù)據(jù)進(jìn)行標(biāo)準(zhǔn)化[31-32]。這種數(shù)據(jù)標(biāo)準(zhǔn)化方法經(jīng)過(guò)處理后數(shù)據(jù)符合標(biāo)準(zhǔn)正態(tài)分布,即均值為0,標(biāo)準(zhǔn)差為1,轉(zhuǎn)化函數(shù)為:
(23)
式中:x為輸入數(shù)據(jù);μ為均值;σ為方差或標(biāo)準(zhǔn)差。深度學(xué)習(xí)模型定義采用Keras框架的序貫?zāi)P?Sequential),首先在模型中添加一個(gè)全連接層(Dense),然后添加5個(gè)隱藏層,最后一層是輸出層,神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)為:3-200-200-100-100-100-1,如圖2所示,即輸入層為3個(gè)神經(jīng)元,對(duì)應(yīng)變量X、Y和深度(Depth),5個(gè)隱藏層分別有200,200,100,100,100個(gè)神經(jīng)元,加入隱藏層可以擬合更加復(fù)雜模型,輸出層為1個(gè)神經(jīng)元,對(duì)應(yīng)預(yù)測(cè)的Q值,這是典型的多元回歸網(wǎng)絡(luò)模型設(shè)置,用來(lái)預(yù)測(cè)單個(gè)連續(xù)值。激活函數(shù)則采用線性整流函數(shù)指代數(shù)學(xué)中的斜坡函數(shù),即Relu函數(shù):f(x)=max(0,x),如圖3所示,加入激活函數(shù)來(lái)擬合非線性模型。模型定義后,則需要對(duì)模型進(jìn)行訓(xùn)練和學(xué)習(xí)。模型學(xué)習(xí)過(guò)程中,最優(yōu)化算子選用自適應(yīng)矩估計(jì)(adaptive moment estimation,ADAM)算子,損失函數(shù)選用最小平方誤差函數(shù)(mean squared error,MSE)。模型訓(xùn)練達(dá)到精度要求后,就形成了一個(gè)三維多元非線性回歸模型,只要輸入工區(qū)內(nèi)任意點(diǎn)的X、Y和深度,就可以預(yù)測(cè)該點(diǎn)的Q值,從而實(shí)現(xiàn)三維Q值建模。算法流程如圖4所示。
圖2 Keras貫序深度學(xué)習(xí)模型結(jié)構(gòu)
圖3 Relu激活函數(shù)
圖4 算法流程
為了便于開(kāi)展微測(cè)井正演分析,結(jié)合生產(chǎn)實(shí)際,建立了如圖5所示的4層水平層狀介質(zhì)模型。正演模擬計(jì)算采用Futteman方程頻率域Q衰減算法實(shí)現(xiàn)。初至波采用非零相位雷克子波,相位為30°,主頻為40Hz。依次在圖5中的1,2,3三個(gè)位置激發(fā),地面接收。為使模擬更加符合實(shí)際情況,在初至后面添加一個(gè)主頻為60Hz,相位為60°的子波,和初至子波呈半疊置關(guān)系,同時(shí)在初至地震信號(hào)添加了隨機(jī)噪聲。
圖5 微測(cè)井正演模型示意
圖6a中藍(lán)色波形是圖5中炮點(diǎn)3激發(fā)地面接收得到的地震記錄,紅色波形是利用CFMP方法匹配得到的初至波形。圖6b是采用CFMP方法提取的初至信號(hào)(紅色)和在炮點(diǎn)3直接激發(fā)不考慮后續(xù)波形和噪聲獲得的初至波(藍(lán)色)疊合顯示。圖6c中藍(lán)色波形是圖5中在炮點(diǎn)2激發(fā)地面接收得到的地震記錄,紅色波形是利用CFMP方法匹配得到的初至波形。圖6d是采用CFMP方法提取的初至信號(hào)(紅色)和在炮點(diǎn)2直接激發(fā)不考慮后續(xù)波形和噪聲獲得的初至波(藍(lán)色)疊合顯示??梢钥吹?盡管有噪聲和續(xù)至波形的影響,采用CFMP方法提取的初至波和實(shí)際初至波基本一致,表明該算法在初至波提取方面具有較好的抗干擾能力。
圖7a中的曲線是在炮點(diǎn)2(藍(lán)色)和炮點(diǎn)3(紅色)激發(fā)地震波地面接收的初至波通過(guò)加窗截取一段信號(hào)后初至信號(hào)求取的頻譜,受噪聲和續(xù)至波的影響,頻譜呈毛刺狀變化。圖7b是用這兩炮的頻譜計(jì)算的譜比值曲線,可以看到,該曲線抖動(dòng)嚴(yán)重。圖7c中藍(lán)虛線是譜比值曲線取對(duì)數(shù)得到的對(duì)數(shù)譜比曲線,同樣曲線抖動(dòng)嚴(yán)重;紅實(shí)線是基于最小二乘線性擬合的結(jié)果,進(jìn)一步利用擬合斜率,由(16)式可以求得第2層的Q值。
圖6 在炮點(diǎn)3激發(fā)合成的微測(cè)井地震記錄(含噪)(a)、采用CFMP方法提取的初至信號(hào)和無(wú)噪聲衰減初至波疊合顯示(b)、在炮點(diǎn)2激發(fā)合成微測(cè)井地震記錄(含噪)(c)以及采用CFMP方法提取的初至信號(hào)和無(wú)噪聲衰減初至波疊合顯示(d)
圖7 常規(guī)的初至波加窗頻譜分析結(jié)果(a)、譜比值曲線(b)以及對(duì)數(shù)譜比曲線與基于最小二乘線性擬合結(jié)果(c)
圖8a中的曲線是在炮點(diǎn)2(藍(lán)色)和炮點(diǎn)3(紅色)激發(fā)地面接收的地震波通過(guò)CFMP提取初至信號(hào)求取的頻譜,由于抗噪性好,計(jì)算的頻譜較為光滑。圖8b是用這兩炮的頻譜計(jì)算的譜比值曲線,可以看到,該曲線基本呈線性。圖8c中藍(lán)虛線是譜比值曲線取對(duì)數(shù)得到的對(duì)數(shù)譜比曲線,紅實(shí)線是基于最小二乘線性擬合的結(jié)果,基于(16)式可以求得第2層的Q值。對(duì)比分析可以看出,與傳統(tǒng)的加窗譜分析方法相比,匹配追蹤得到的頻譜能夠分離干擾信號(hào),得到的頻譜更加平滑,頻譜分析精度更高,同時(shí)加入正則化的譜比相較于傳統(tǒng)的譜比法運(yùn)算更加穩(wěn)定。
表1給出了模型Q值、估計(jì)Q值及其誤差。從表1可以看出,CFMP譜比法和常規(guī)譜比法相比,抗噪能力強(qiáng),Q值估算累計(jì)誤差小,精度高。
圖8 采用CFMP方法提取的初至信號(hào)頻譜分析結(jié)果(a)、譜比值曲線(b)以及對(duì)數(shù)譜比曲線與基于最小二乘線性擬合結(jié)果(c)
表1 模型Q值、估計(jì)Q值及其誤差
實(shí)測(cè)微測(cè)井?dāng)?shù)據(jù)來(lái)自中國(guó)石油華北油田。探區(qū)地表平緩,但表層結(jié)構(gòu)松散且復(fù)雜多變、速度低,對(duì)地震波的吸收、衰減嚴(yán)重,地震資料分辨率低。巖石物理分析發(fā)現(xiàn),介質(zhì)的吸收作用與介質(zhì)孔隙度、固結(jié)程度以及孔隙充填物關(guān)系密切,所以表層未固結(jié)巖層的吸收規(guī)律與深層固結(jié)巖層的規(guī)律不同。圖9為由近地表取樣實(shí)驗(yàn)室測(cè)定的泥、黃土和沙3種巖樣Q值與速度交會(huì)分析結(jié)果。由于表層通常是砂、泥混合的,可以看出,近地表的Q值為2~15,速度為300~800m/s。
圖9 近地表不同巖石樣本測(cè)試Q值與速度交會(huì)分析結(jié)果
圖10給出了探區(qū)內(nèi)微測(cè)井測(cè)點(diǎn)位置,探區(qū)面積約60km2,區(qū)內(nèi)共有42個(gè)微測(cè)井測(cè)量點(diǎn)。圖11是采用CFMP對(duì)數(shù)譜比法估算的圖10中不同點(diǎn)所在位置的Q值函數(shù)。
圖10 工區(qū)微測(cè)井測(cè)點(diǎn)分布
圖11 采用CFMP對(duì)數(shù)譜比法估算的圖10中不同點(diǎn)所在位置的Q值函數(shù)
以探區(qū)內(nèi)42口微測(cè)井資料估計(jì)的近地表Q值函數(shù)中的36口井的數(shù)據(jù)構(gòu)成訓(xùn)練集作為深度學(xué)習(xí)模型的輸入樣本,即訓(xùn)練標(biāo)簽,其余6口井作為測(cè)試集。進(jìn)一步按照前述深度學(xué)習(xí)三維Q值建模原理,通過(guò)深度神經(jīng)網(wǎng)絡(luò)的多元非線性回歸,可得到一個(gè)三維近地表Q值回歸數(shù)學(xué)模型。圖12給出了圖10中6號(hào)井點(diǎn)微測(cè)井資料估算的Q值函數(shù)(測(cè)試集中的一個(gè)樣本)與深度學(xué)習(xí)Q值回歸數(shù)學(xué)模型所得預(yù)測(cè)值的對(duì)比結(jié)果,可以看出,預(yù)測(cè)值與微測(cè)井估算Q值之間的絕對(duì)誤差值均介于±0.1內(nèi),基本可以忽略,說(shuō)明二者具有良好的一致性。
圖13為基于深度學(xué)習(xí)多元非線性回歸所得近地表Q值數(shù)學(xué)模型建立的研究區(qū)三維近地表Q值模型。圖14a和圖14b分別為10m和20m的近地表Q值等深度切片。由圖14可見(jiàn),Q值空間變化規(guī)律性較強(qiáng),盡管控制樣本點(diǎn)很少,但是訓(xùn)練獲得的多元非線性回歸模型預(yù)測(cè)的Q值變化較為平滑,證明了該算法在數(shù)據(jù)量較為稀疏的條件下可以很好地構(gòu)建三維模型。圖14中紅線是三維工區(qū)的邊界,邊界線外預(yù)測(cè)的結(jié)果可靠性要差一些。
圖12 6號(hào)井點(diǎn)微測(cè)井資料估算的Q值函數(shù)(a)、深度學(xué)習(xí)Q值回歸數(shù)學(xué)模型所得的預(yù)測(cè)值(b)和預(yù)測(cè)誤差(c)
圖13 近地表三維Q值模型
圖14 近地表Q值深度切片
地層非彈性吸收引起的地震波能量衰減的量可通過(guò)反Q濾波方法對(duì)衰減量進(jìn)行補(bǔ)償校正,以提高地震剖面的分辨率。在反Q濾波的過(guò)程中,外推步長(zhǎng)是地震道的時(shí)間采樣率,并要求在每一個(gè)采樣點(diǎn)上輸出結(jié)果,這樣計(jì)算量巨大,限制了其在疊前處理過(guò)程中的應(yīng)用。根據(jù)地震波在地下介質(zhì)中傳播的衰減規(guī)律,可將近地表層Q值轉(zhuǎn)化為等效Q值:
(24)
式中:Qj為近地表第j層的Q值;Q0N(τ)表示從地表到第N層的等效Q值,描述了地面激發(fā)點(diǎn)與地下相應(yīng)點(diǎn)間的Q值關(guān)系。最后將單程旅行時(shí)概念下的等效Q值轉(zhuǎn)化為雙程旅行時(shí)概念下的等效Q(τ)值?;诮乇淼刃(τ)值的衰減公式表示如下:
(25)
式中:μ(τ)=1/πQ(τ);A(0,ω)是地面接收的地震信號(hào)振幅;A(τ,ω)是傳過(guò)等效Q層后的地震信號(hào)振幅;ω是圓頻率;ωr是參考圓頻率;i是虛數(shù)單位。該方法可以對(duì)振幅和相位同時(shí)進(jìn)行校正,且等效Q值的引入更好地適應(yīng)了近地表的衰減補(bǔ)償處理,更有利于提高地震資料分辨率。按上述方法可以將隨深度變化的近地表層Q模型分為多個(gè)等效Q模型,在每一等效Q層內(nèi)實(shí)施反Q濾波。將波場(chǎng)從地面延拓到第N層的頂部,可以遞歸實(shí)現(xiàn),即:
A(τ0,ω)→A(τ1,ω)→…→A(τn-1,ω)
(26)
在近地表反Q濾波中,則是將地面地震記錄按反Q濾波的形式延拓到近地表Q值模型的底部,在近地表模型底部以下,給定一個(gè)較大的Q值(例如200),以保證反Q濾波只受近地表參數(shù)的影響。
圖15為近地表Q值補(bǔ)償濾波前、后的單炮記錄,顯然,近地表Q值補(bǔ)償后地震記錄的分辨率有明顯的提高。圖16是圖15所示的單炮記錄近地表Q值補(bǔ)償前、后的頻譜,可以看出,無(wú)論是1500ms以上還是以下,高頻成分都得到了有效補(bǔ)償,頻帶展寬明顯。圖17為近地表Q值補(bǔ)償前、后的疊前時(shí)間偏移剖面,其中,紅色波形是插入的合成地震記錄,可以看出,近地表Q值補(bǔ)償后,目標(biāo)層分辨率明顯提高,與合成地震記錄的一致性更好。
圖15 近地表Q值補(bǔ)償前(a)、后(b)的單炮記錄
圖16 圖15所示的單炮記錄近地表Q值補(bǔ)償前、后的頻譜
圖17 近地表Q值補(bǔ)償前(a)、后(b)的疊前時(shí)間偏移剖面
本文提出了一種基于非零相位雷克子波的復(fù)數(shù)域快速匹配追蹤分解并結(jié)合對(duì)數(shù)譜比法估算微測(cè)井?dāng)?shù)據(jù)近地表Q值的方法,并基于深度學(xué)習(xí)構(gòu)建了多元非線性回歸算子,實(shí)現(xiàn)了近地表Q值三維建模,進(jìn)而將Q值模型用于實(shí)際資料處理,得出以下結(jié)論。
1) 采用基于非零相位雷克子波原子庫(kù)改進(jìn)的快速?gòu)?fù)數(shù)域匹配追蹤方法,可自動(dòng)提取微測(cè)井初至波的頻譜信息和初至?xí)r間,提高了初至波頻譜求取的精度和抗噪能力,并提高了微測(cè)井資料處理的工作效率。
2) 通過(guò)引入整形正則化算子和最小平方反演方法求取譜比值,提高了譜比算法的準(zhǔn)確性和穩(wěn)定性。
3) 將工區(qū)內(nèi)不同位置的微測(cè)井?dāng)?shù)據(jù)估計(jì)得到的Q值函數(shù)作為已知樣本標(biāo)簽,基于遷移學(xué)習(xí)原理,通過(guò)深度學(xué)習(xí)訓(xùn)練形成近地表Q值的多元非線性回歸算子,可建立較傳統(tǒng)方法更精確的三維近地表Q值模型。
4) 華北油田實(shí)際資料應(yīng)用結(jié)果表明,本文近地表Q值估算方法具有合理性和準(zhǔn)確性,且高效、抗噪能力強(qiáng)。利用建立的三維近地表Q值模型進(jìn)行近地表Q值補(bǔ)償,可展寬目的層反射波頻帶,提高地震資料分辨率。
研究發(fā)現(xiàn)盡管可以通過(guò)匹配追蹤獲取精確的初至波頻譜和旅行時(shí),但是在譜比計(jì)算中,需要選取計(jì)算譜比值的頻帶寬度參數(shù),而估算的Q值對(duì)譜比值比較敏感,因此往往會(huì)造成估算結(jié)果的不穩(wěn)定??梢圆捎美碚摼认鄬?duì)差一點(diǎn)的質(zhì)心頻移法代替譜比法進(jìn)行Q值估算,以提高算法的穩(wěn)定性。