本發(fā)明屬于轉(zhuǎn)換波地震資料處理技術(shù)領(lǐng)域,具體來講,涉及一種高分辨率轉(zhuǎn)換波裂縫預(yù)測方法。
背景技術(shù):
目前,相干體裂縫特征描述是探尋裂縫性油氣藏的一項重要手段,而獲取裂縫預(yù)測所需高分辨率數(shù)據(jù)體、裂縫參數(shù)的求取一直是裂縫油氣藏勘探中重點研究的問題,尤其轉(zhuǎn)換波低信噪比和分辨率的特點,需要提高其分辨率,便于后續(xù)解釋的薄層分辨。但是,最小相位和白噪假設(shè)為前提的常規(guī)高分辨率處理方法,在主頻帶范圍限制下,提高轉(zhuǎn)換波分辨率的能力有限,無法達到后續(xù)裂縫預(yù)測分辨多薄層最小厚度的要求,所以需要采取新的高分辨率處理思路,即:轉(zhuǎn)換波高精度反射系數(shù)求取與高頻子波高分辨率重構(gòu)的新思路,該思路不受任何假設(shè)和頻帶限制,能完成常規(guī)手段無法分辨多薄層最小厚度的要求,獲取裂縫預(yù)測滿意的高品質(zhì)剖面。
通常,反射系數(shù)的求取需要構(gòu)建多薄層目標函數(shù),然而,由于常規(guī)方法構(gòu)建的目標函數(shù)存在諸多不合理之處,故而,在求取過程中的反射系數(shù)存在極值和位置、大小上的誤差,這種誤差導(dǎo)致后續(xù)高分辨率剖面的錯像。此外,目標函數(shù)的求解方法有很多,比如模擬退火法、共軛梯度法、匹配追蹤法等這些算法,雖然都能求取目標函數(shù)的反射系數(shù),但精度和穩(wěn)定性較差,且對新構(gòu)建的目標函數(shù)適應(yīng)性較差。
三維相干體技術(shù)是20世紀90年代后期興起的一項十分有效的地震解釋技術(shù),該技術(shù)主要從相鄰地震道相互之間的相干性出發(fā),給出一種定量描述,它能突出斷層處地震波的變化情況,形成相干值異常區(qū)域,從而預(yù)測斷裂展布。M.Bahorich等提出在振幅上進行互相關(guān)相干分析算法,但該算法計算出的斷點非相干性點模糊,噪聲重、能量弱。K.J.Marfurt年提出基于相似系數(shù)相干算法具有穩(wěn)定性好、斷層刻畫精度高的特點,但橫向分辨率低的劣勢無法改變。Gersztenkorn等提出基于本征結(jié)構(gòu)分析的算法,利用Hilbert變換計算相干體提高其穩(wěn)定性和抗噪性,但該算法用于高陡復(fù)雜地區(qū)時,如不考慮傾角則不如第二代算法效果好,如考慮傾角則計算效率無法應(yīng)用于生產(chǎn)中。
技術(shù)實現(xiàn)要素:
本發(fā)明的目的在于解決現(xiàn)有技術(shù)存在的上述不足中的至少一項。例如,本發(fā)明的目的之一在于提供一種具有更高分辨率的轉(zhuǎn)換波裂縫預(yù)測方法,以使得提取的裂縫預(yù)測剖面分辨率更高、裂縫走向更加清晰、細節(jié)刻畫更加明顯。
為了實現(xiàn)上述目的,本發(fā)明提供一種高分辨率轉(zhuǎn)換波裂縫預(yù)測方法。所述方法包括:采集轉(zhuǎn)換波地震數(shù)據(jù),進行預(yù)處理和疊前偏移后得到轉(zhuǎn)換波偏移剖面,在轉(zhuǎn)換波偏移剖面上取n道地震記錄x(x,t),x=1,2,…,n;在偏移剖面上建立稀疏約束和波阻抗約束的多薄層目標函數(shù);對多薄層目標函數(shù)利用頻域映射與L1范數(shù)聯(lián)合優(yōu)化算法求取高精度反射系數(shù);構(gòu)建高頻子波同時保持該高頻子波具有確定的頻率寬度,用高頻子波與高精度反射系數(shù)得到高分辨率轉(zhuǎn)換波偏移剖面;利用改進的第三代特征值相干體計算法,對高分辨率轉(zhuǎn)換波偏移剖面進行相干數(shù)據(jù)體切片提取,從而實現(xiàn)裂縫預(yù)測;其中,
所述目標函數(shù)為下式(1):
其中,r表示反射系數(shù),ae代表反射系數(shù)的偶分量比例,ao代表反射系數(shù)的奇分量比例,re為反射系數(shù)的偶分量,ro為反射系數(shù)的奇分量,A為變換矩陣,b為地震數(shù)據(jù)和子波相關(guān)矩陣,λ為稀疏因子,為阻抗因子,C為積分算子矩陣,ξp為縱波波阻抗矩陣。
與現(xiàn)有技術(shù)相比,本發(fā)明的有益效果包括:能夠基于轉(zhuǎn)換波地震數(shù)據(jù)構(gòu)建更加合理的多薄層目標函數(shù),并能有效提高目標函數(shù)求解精度和抗噪性,能夠獲得高精度反射系數(shù),能夠使得提取的裂縫預(yù)測剖面分辨率更高、裂縫走向更加清晰、細節(jié)刻畫更加明顯。
附圖說明
圖1示出了通過現(xiàn)有的預(yù)處理和疊前偏移后得到的轉(zhuǎn)換波偏移剖面。
圖2示出了根據(jù)本發(fā)明的方法求取高精度反射系數(shù)對應(yīng)的剖面。
圖3示出了根據(jù)本發(fā)明的方法得到的高分辨率轉(zhuǎn)換波偏移剖面。
圖4A示出了圖1的轉(zhuǎn)換波偏移剖面的振幅譜數(shù)據(jù);圖4B示出了根據(jù)本發(fā)明的方法得到的高分辨率轉(zhuǎn)換波偏移剖面的振幅譜數(shù)據(jù)。
圖5A示出了利用改進的第三代特征值相干體計算法對圖1的轉(zhuǎn)換波偏移剖面進行相干切片提取的效果圖;圖5B示出了根據(jù)本發(fā)明的方法對高分辨率轉(zhuǎn)換波偏移剖面進行相干切片提取的效果圖。
具體實施方式
在下文中,將結(jié)合附圖和示例性實施例來詳細說明本發(fā)明的高分辨率轉(zhuǎn)換波裂縫預(yù)測方法。
在本發(fā)明的一個示例性實施例中,高分辨率轉(zhuǎn)換波裂縫預(yù)測方法可通過下面的步驟(1)至步驟(5)來實現(xiàn):
(1)對采集到的轉(zhuǎn)換波地震數(shù)據(jù),經(jīng)過常規(guī)預(yù)處理和疊前偏移后得到的轉(zhuǎn)換波偏移剖面,取n道地震記錄x(x,t),x=1,2,…,n;其中,x表示道數(shù),t表示時間。
(2)在偏移剖面上建立由稀疏約束和波阻抗約束雙重約束的多薄層目標函數(shù);
步驟(2)可按下列方式完成:
(a)對于不含噪聲的褶積模型,地震合成記錄在頻率域中可表示為:
S(f)=W(f)·R(f) (1)
其中,S(f)表示頻率域下地震記錄(即,上面的x(x,t)),W(f)表示頻率域下地震子波,R(f)表示頻率域下反射系數(shù);f為頻率,時間t經(jīng)過傅里葉變換后變換到頻率f中。
利用頻率域下地震記錄與地震子波的豐富信息,利用合適的常規(guī)反演算法(如,共軛梯度法、匹配追蹤法、模擬退火等)消除子波(即,地震子波)影響(如,子波主頻、旁瓣干擾等會影響反演求取反射系數(shù)精度的內(nèi)容),從記錄中得到所需的反射系數(shù)(這里,為了得到反射系數(shù)R(f))。利用反射系數(shù)奇偶分解原理,可得到頻率域下目標函數(shù)表達式:
其中,ae代表反射系數(shù)R(f)的偶分量比例,ao代表反射系數(shù)R(f)的奇分量比例;Re表示R(f)的實部,Im表示R(f)的虛部。
(b)為建立更符合實際地層情況的反射系數(shù)模型,對步驟(a)中的多薄層反射系數(shù)模型下的目標函數(shù)進行推導(dǎo)或改進。
由奇偶反射系數(shù)分解原理,得多薄層反射系數(shù)奇偶分量可表示為:
re(t)代表反射系數(shù)r(t)的偶分量,ro(t)代表反射系數(shù)r(t)的奇分量,t表示時間,T表示時間厚度;r(t)為時間域反射系數(shù),r(t)經(jīng)過傅里葉變換后到上面的頻率域反射系數(shù)R(f);g(t)當前分析位置反射系數(shù),ge(t)代表g(t)的偶分量,go(t)代表g(t)的奇分量。
對式(3)進行傅里葉變換,得到其頻率域?qū)嵅颗c虛部表達式:
其中,π為圓周率。
利用式(2)和式(4)得到多薄層反射模型的目標函數(shù)表達式為:
其中,re代表反射系數(shù)r(t)的偶分量,ro代表反射系數(shù)r(t)的奇分量,Δt表示時移量;i無物理意義,是指數(shù)與正弦函數(shù)變化標識符。
為了便于求解,將式(5)可寫成矩陣形式:
其中,
其中,t為時間,t1、t2、……、tN分別表示第1、2、……、N個時間,f為頻率,f1、f2、……、fM分別表示第1、2、……、M個分析頻率范圍內(nèi)的頻率,T1、T2……、TN/2分別表示第1、2、……、N/2個時間厚度,Δt表示時移量
(c)由于該目標函數(shù)(6)(即,式(6))有多個全局最優(yōu)解,為了降低求解問題時的多解性和不確定性,需要加上稀疏化約束,將目標函數(shù)改為范數(shù)求解形式為:
其中,λ為稀疏因子,λ一般在0-1范圍內(nèi)取值,|| ||1表示L1范數(shù),表示歐幾里得(Euclidean)范數(shù)的平方。
(d)由于波阻抗可看作是反射系數(shù)對時間的積分,因此發(fā)明人把它引入到目標函數(shù)中作為反射系數(shù)r的先驗約束。
由于轉(zhuǎn)換波相對于縱波來說其信噪比和分辨率較低,為了提高反射系數(shù)精度,本發(fā)明使用分辨率更高的縱波波阻抗作為轉(zhuǎn)換波的約束條件。
若I(t)為縱波波阻抗值,則縱波相對波阻抗εt為:
其中,I(t0)表示縱波初始波阻抗值,r為反射系數(shù),η表示在時間t0與時間t中取值符合,t0表示初始時間,tM-1表示結(jié)束時間。
把上式(8)寫成矩陣形式:
ξp=Cr (9)
其中C為積分算子矩陣,其離散形式可表示為:
則由最小平方定義波阻抗約束得到新的目標函數(shù)為:
其中,r表示反射系數(shù),ae和ao分別代表反射系數(shù)r(t)的偶奇分量比例,re是反射系數(shù)的偶分量,ro是反射系數(shù)的奇分量;A為變換矩陣,b是地震數(shù)據(jù)和子波相關(guān)矩陣,λ是稀疏因子,是阻抗因子(其大于零即可),C為積分算子矩陣,ξp是縱波波阻抗矩陣)
(3)對于步驟(2)中的式(11),利用頻域映射與L1范數(shù)聯(lián)合優(yōu)化算法求取高精度反射系數(shù)。
本步驟(3)可按下列方式完成:
(a)對式(11)進行求解中,由于波阻抗約束項是由測井信號得到的已知先驗約束條件,所以式(11)的求解類似于求矩陣方程Ax=b在含有噪聲情況下的稀疏解。
可將目標函數(shù)(11)的求解問題轉(zhuǎn)化為下式:
其中,υ為極小權(quán)重值,G為復(fù)數(shù)集合;本式(12)中,A表示Ae時b表示be,X表示re;A表示Ao時b表示bo,X表示ro。
如果引入變量q∈Gm,那么式(12)等價于式(13):
式(13)對應(yīng)的增廣拉格朗日子問題可表示為:
其中,yT表示拉格朗日乘子的共軛轉(zhuǎn)置運算,β為罰參數(shù)。
利用式(14)進行頻域映射與L1范數(shù)聯(lián)合優(yōu)化算法實現(xiàn)步驟如下:
1)令k=0對qk,xk,yk賦初始值,并給定υ,β常數(shù)值(這里,υ,β為人為給定,可為大于零的常值),然后進行以下步驟2)至5)的算法運算,如果滿足終止準則(例如,終止準則可以為:最優(yōu)解值或者迭代次數(shù)),則算法完成運算;否則,進行步驟2);
2)令x=xk,y=y(tǒng)k,求解r的子問題得到qk+1:
3)令q=qk+1,y=y(tǒng)k則關(guān)于x的極小化問題式(14)等價于式(16):
那么可以通過式(17)近似求解來完成對式(16)的精確求解:
其中,τ為大于零的鄰近參數(shù)(這里,τ為人工給定,其大于零即可),gk為x=xk時二次項的梯度,那么式(17)通過下式(18)求解得到xk+1:
其中,ο表示逐元素相乘,Shrink(,)表示一維收縮算子,sign表示符號函數(shù)。
4)令x=xk+1,q=qk+1,求解y的子問題得到y(tǒng)k+1:
yk+1=y(tǒng)k-γβ(Axk+1`+qk+1-b) (19)
其中,γ為大于零的常數(shù)(這里,γ人工給定,其大于零即可)。
5)令k=k+1,重復(fù)步驟1)~步驟4)進行迭代運算。
(4)基于步驟(3)反射系數(shù)與高頻子波重構(gòu),得到高分辨率轉(zhuǎn)換波偏移剖面。
構(gòu)建高頻子波的同時需要保持該高頻子波具有一定的頻帶寬度,以使得高分辨率轉(zhuǎn)換波偏移剖面是寬頻的??衫贸R?guī)的精度較高的俞式子波建立高頻子波信號,亦可用諸如雷克子波等零相位子波建立高頻子波信號。
時間域中俞氏子波表達式為:
其中q、p分別表示子波的低頻值和高頻值,t為時間。
用式(20)的子波與步驟(3)求取的高精度反射系數(shù)得到高分辨率轉(zhuǎn)換波偏移剖面。
(5)對步驟(4)的高分辨率轉(zhuǎn)換波偏移剖面進行改進的第三代特征值相干算法C3相干體切片裂縫預(yù)測。
相干數(shù)據(jù)體技術(shù)主要利用相鄰地震道的相似性原理,描述地層及巖性的橫向不均勻性。通過相干數(shù)據(jù)體算法,對地震數(shù)據(jù)體的不連續(xù)性進行系統(tǒng)分析,使解釋人員能有效識別構(gòu)造和斷層的分布情況,避免解釋的隨意性。
(a)為了盡可能壓制噪聲干擾,提高相干體計算數(shù)據(jù)的分辨率,Gersztenkorn等提出了第三代特征值相干體計算法。
首先,利用Marfurt等在第二代多道相似性相干體計算中定義的橢圓和矩形分析窗口中,以坐標為(x,y)的分析點為中心,時窗中心t=nΔt,則時間窗口內(nèi)含2k+1個元素的分析道與窗口內(nèi)的J道數(shù)據(jù)排列成矩陣:
其中,表示第j道在時間t=mΔt-pΔxj-qyj沿著視傾角(p,q)插值后的振幅值。
那么矩陣U的第m列是分析體內(nèi)每個地震道在時間t=mΔt-pΔxj-qyj(j=1,2,…,J)的振幅值,則沿著一對視傾角(p,q),時窗中心在t=mΔt的分析窗口內(nèi)的2k+1個采樣點的數(shù)據(jù)協(xié)方差矩陣為:
如λj為分析窗口內(nèi)地震數(shù)據(jù)協(xié)方差矩陣式(22)的特征值,且有λ1≥λ2≥…≥λJ。
則在以分析點為中心的分析窗口內(nèi),沿著視傾角(p,q)的J道c3特征值相干計算算法式可表示為:
則在以分析點為中心的分析窗口內(nèi),沿視傾角(p,q=0)的J道c3特征值相干計算算法式可表示為:
則沿著所有傾角的c3特征值相干計算算法式可表示為:
(b)因為第三代特征值相干算法C3是在子空間中計算,比前兩代相干算法更能很好地壓制噪聲,但其計算量大、效率慢的缺點卻影響了其在實際生產(chǎn)中的經(jīng)濟效益。
對該算法進行常規(guī)改進。例如,可利用基于相似性算法的高效率性估算各點的視傾角和方位角,然后在大于10倍的分析窗口內(nèi)計算視傾角和方位角的平均值再在以分析點為中心的分析窗口內(nèi),沿著區(qū)域平滑后的視傾角計算地震數(shù)據(jù)協(xié)方差矩陣的特征值及相干值:
利用式(26)對步驟(4)的高分辨率剖面進行相干數(shù)據(jù)體切片提取,從而進行裂縫預(yù)測。
圖1示出了通過現(xiàn)有的預(yù)處理和疊前偏移后得到的轉(zhuǎn)換波偏移剖面,其橫坐標表示CMP道號,圖1中也示出了測井系數(shù)(即,圖1中位于道號791和841之間的縱向曲線)。圖2示出了根據(jù)本發(fā)明的方法求取高精度反射系數(shù)對應(yīng)的剖面。由圖2可以看出,本發(fā)明的方法求取高精度反射系數(shù)與測井系數(shù)(縱向曲線)匹配較好,各層位置和極性求取較佳。
圖3示出了根據(jù)本發(fā)明的方法得到的高分辨率轉(zhuǎn)換波偏移剖面??梢钥闯?,圖3明顯比圖1的分辨率高,且淺、中、深層的同相軸(例如,800ms附近的淺層、1000ms附近的中層、1150ms附近的深層的同相軸)連續(xù)性更好、斷點更清晰、斷層更明了,便于后續(xù)解釋人員的裂縫預(yù)測處理。
圖4A示出了圖1的轉(zhuǎn)換波偏移剖面的振幅譜數(shù)據(jù);圖4B示出了根據(jù)本發(fā)明的方法得到的高分辨率轉(zhuǎn)換波偏移剖面的振幅譜數(shù)據(jù)。由圖4A和圖4B可知,經(jīng)本發(fā)明的方法處理而得到的高分辨率轉(zhuǎn)換波偏移剖面的振幅譜數(shù)據(jù)的頻帶得到明顯拓寬,而且其振幅譜低頻信息得到完整的保留。
圖5A示出了利用改進的第三代特征值相干體計算法對圖1的轉(zhuǎn)換波偏移剖面進行相干切片提取的效果圖;圖5B示出了根據(jù)本發(fā)明的方法對高分辨率轉(zhuǎn)換波偏移剖面進行相干切片提取的效果圖。由圖5A和圖5B可知,經(jīng)本發(fā)明的方法處理后的切片相干效果更好,裂縫刻畫細節(jié)更清晰,避免了低分辨率數(shù)據(jù)為后續(xù)解釋人員帶來的假象。
綜上所述,本發(fā)明通過利用稀疏約束和波阻抗雙重約束構(gòu)建多薄層目標函數(shù),提高目標函數(shù)求解精度和抗噪性;利用頻域映射與L1范數(shù)聯(lián)合優(yōu)化算法求取高精度反射系數(shù),并通過與高頻子波重構(gòu)得到分辨率更高的新地震剖面,改善轉(zhuǎn)換波的分辨率效果,利用改進的C3相干體切片進行裂縫預(yù)測,可有效提高地震勘探裂預(yù)測和描述地下構(gòu)造的能力。也就是說,本發(fā)明的高分辨率轉(zhuǎn)換波裂縫預(yù)測方法能夠使得提取的裂縫預(yù)測剖面分辨率更高、裂縫走向更加清晰、細節(jié)刻畫更加明顯。
盡管上面已經(jīng)結(jié)合附圖和示例性實施例描述了本發(fā)明,但是本領(lǐng)域普通技術(shù)人員應(yīng)該清楚,在不脫離權(quán)利要求的精神和范圍的情況下,可以對上述實施例進行各種修改。