陳可洋,吳沛熹,楊 微
(1.中國(guó)石油大慶油田有限責(zé)任公司勘探開(kāi)發(fā)研究院,黑龍江大慶163712;2.中國(guó)石油大慶油田有限責(zé)任公司第六采油廠(chǎng),黑龍江大慶163712)
擴(kuò)散濾波方法在地震資料處理中的應(yīng)用研究
陳可洋1,吳沛熹2,楊 微2
(1.中國(guó)石油大慶油田有限責(zé)任公司勘探開(kāi)發(fā)研究院,黑龍江大慶163712;2.中國(guó)石油大慶油田有限責(zé)任公司第六采油廠(chǎng),黑龍江大慶163712)
從各向異性擴(kuò)散濾波方程出發(fā),給出了其對(duì)應(yīng)的數(shù)值離散公式及其最簡(jiǎn)單的各向同性擴(kuò)散濾波離散計(jì)算公式,開(kāi)展了擴(kuò)散濾波方法在低頻逆時(shí)偏移噪音的多尺度分離、地震反射波和散射波波場(chǎng)的相對(duì)分離、保斷層邊緣的隨機(jī)噪音壓制以及地震數(shù)據(jù)插值處理等4個(gè)方面的應(yīng)用,并給出了其對(duì)應(yīng)的應(yīng)用假設(shè)條件。計(jì)算結(jié)果表明,擴(kuò)散濾波方法可在不同應(yīng)用假設(shè)條件下取得較好的處理效果。該方法可指導(dǎo)實(shí)際地震資料數(shù)字處理。
擴(kuò)散濾波;噪音壓制;數(shù)據(jù)插值;波場(chǎng)分離;地震勘探
基于偏微分方程的信號(hào)處理技術(shù)于20世紀(jì)90年代初興起,該技術(shù)能夠保留圖像的細(xì)節(jié)特征,并檢測(cè)圖像的局部特征,能夠避免傳統(tǒng)的圖像處理方法(高斯濾波和中值濾波等)在去噪時(shí)會(huì)模糊甚至破壞圖像細(xì)節(jié)信息的缺陷,因此,它在圖像處理和地震勘探等多個(gè)工程領(lǐng)域中得到了廣泛應(yīng)用[1-2]。其中最具代表性的是擴(kuò)散濾波方法,其機(jī)理來(lái)源于熱的擴(kuò)散現(xiàn)象,它將輸入待處理的圖像數(shù)據(jù)作為初始的邊界條件,通過(guò)求解擴(kuò)散濾波方程得到預(yù)定擴(kuò)散時(shí)間后的圖像。在擴(kuò)散濾波迭代過(guò)程中,可根據(jù)圖像的局部特征施加一定的約束條件,從而實(shí)現(xiàn)可控制的擴(kuò)散濾波處理。Perona等[3]首次提出用于熱力學(xué)傳導(dǎo)的偏微分方程;Fehmers等[4]將該方法引入到地震資料的處理和解釋中,提出構(gòu)造約束與邊緣保持的各向異性擴(kuò)散濾波方法,用于增強(qiáng)地震數(shù)據(jù)的結(jié)構(gòu)特征,使其易于解釋?zhuān)籐avialle等[5]提出斷層保持的擴(kuò)散濾波方法,使地震數(shù)據(jù)中的斷層更為突出;Kadlec等[6]提出保持層序特征的各向異性擴(kuò)散方法,用于增強(qiáng)地震數(shù)據(jù)中河道砂體特征的連續(xù)性,提高河道自動(dòng)識(shí)別和分割能力。孫夕平等[7]和王緒松等[8]探索了非線(xiàn)性各向異性擴(kuò)散中的一致性增強(qiáng)擴(kuò)散濾波技術(shù)在二維地震剖面保邊濾波中的應(yīng)用。張爾華等[9]對(duì)三維地震數(shù)據(jù)非線(xiàn)性各向異性擴(kuò)散濾波方法進(jìn)行了初步探討,通過(guò)控制斷層和河道等不連續(xù)性區(qū)域的濾波程度,使各向異性擴(kuò)散濾波器具有較好的保邊處理性能。陳可洋等[10]提出了基于擴(kuò)散濾波的地震反射波和散射波的相對(duì)波場(chǎng)分離方法,并取得了較好的應(yīng)用效果。
筆者在前人研究的基礎(chǔ)上,從各向異性擴(kuò)散濾波偏微分方程出發(fā),推導(dǎo)各向異性和各向同性擴(kuò)散濾波的數(shù)值離散公式,總結(jié)當(dāng)前擴(kuò)散濾波方法在低頻逆時(shí)偏移噪音壓制、反射波和散射波波場(chǎng)分離、隨機(jī)噪音壓制和地震數(shù)據(jù)插值處理等4個(gè)重要方面的應(yīng)用情況及假設(shè)條件,以期為實(shí)際地震資料處理和解釋等提供指導(dǎo)。
Perona等[3]提出了用于圖像增強(qiáng)的偏微分方程,即擴(kuò)散濾波方程
式中:t為擴(kuò)散時(shí)間(可不考慮時(shí)間步長(zhǎng)dt,常用迭代次數(shù)N代替),s;div為散度算子,無(wú)量綱;為梯度算子,無(wú)量綱;U為t時(shí)刻的擴(kuò)散濾波結(jié)果,其中U0為t=0時(shí)刻的原始數(shù)據(jù),即擴(kuò)散濾波迭代計(jì)算的初始條件,無(wú)量綱;g(·)為擴(kuò)散函數(shù),Perona等[3]建議其形式為和無(wú)量綱。
以二維地震數(shù)據(jù)為例,在3×3網(wǎng)格窗口內(nèi)對(duì)式(1)進(jìn)行離散化[11],設(shè)窗口中心像素為U(x+2,z+2),則有
式中:gi,k為擴(kuò)散函數(shù)g(·)的二維離散化形式。
對(duì)式(2)進(jìn)行整理可得
令
則有
式(5)為二維各向異性擴(kuò)散濾波方程在3×3網(wǎng)格窗口的離散化形式。此時(shí)的離散化形式分別為
式中Δ:G為擴(kuò)散系數(shù)的歸一化因子,無(wú)量綱;K為梯度U的濾波門(mén)檻,無(wú)量綱。
當(dāng)擴(kuò)散函數(shù)為某一常數(shù)g0時(shí),有??Ut=g0Δ2U,此時(shí)式(5)的離散方程形式為
2.1 低頻逆時(shí)偏移噪音的有效壓制
應(yīng)用假設(shè)條件:低頻逆時(shí)噪音是由炮點(diǎn)波場(chǎng)和檢波點(diǎn)波場(chǎng)中不相關(guān)的波沿傳播路徑作相關(guān)運(yùn)算形成,最終在逆時(shí)成像剖面上疊加的較強(qiáng)能量的背景噪音[12-13]。
對(duì)逆時(shí)偏移剖面分析認(rèn)為,該種噪音具有頻率低和全局平滑的特點(diǎn),而有效的地層信號(hào)掩蓋于此平滑的背景噪音中。擴(kuò)散濾波方法可通過(guò)調(diào)節(jié)擴(kuò)散系數(shù)和迭代次數(shù)實(shí)現(xiàn)相對(duì)低頻的背景噪音與有效地層信號(hào)的分離。
圖1(a)為某地區(qū)實(shí)際地震資料的疊前逆時(shí)偏移疊加剖面,采用的是相關(guān)法逆時(shí)成像。分析圖1(a)可知,低頻噪音較強(qiáng),掩蓋了有效的地層細(xì)節(jié),不利于地層構(gòu)造等特征的有效刻畫(huà),但在對(duì)應(yīng)的深度域波數(shù)譜中,低波數(shù)端能量不突出。圖1(b)為圖1(a)經(jīng)擴(kuò)散濾波方法處理后的逆時(shí)偏移剖面,可見(jiàn)低頻噪音能量得到了有效壓制,地層細(xì)節(jié)特征刻畫(huà)得更加清晰。在其深度域波數(shù)譜中,低波數(shù)端能量有一定的降低,同時(shí)高波數(shù)端能量保持不變。圖1(c)為經(jīng)擴(kuò)散濾波方法處理后被壓制掉的低頻噪音剖面,反映的是圖1(a)背景能量的變化,且在其深度域波數(shù)譜中,低頻能量?jī)H出現(xiàn)于低波數(shù)端。由此可見(jiàn),擴(kuò)散濾波方法可實(shí)現(xiàn)低頻逆時(shí)偏移噪音的有效壓制,從而達(dá)到提高刻畫(huà)地層細(xì)節(jié)的能力。同時(shí),通過(guò)調(diào)節(jié)擴(kuò)散系數(shù)和迭代次數(shù)還可獲得不同尺度條件下的低頻噪音壓制結(jié)果及殘差剖面。
圖1 擴(kuò)散濾波方法在低頻逆時(shí)偏移噪音壓制中的應(yīng)用Fig.1 App lication of diffusion filteringmethod to low-frequency reverse-tim em igration noisesuppression
2.2地震反射波和散射波波場(chǎng)相對(duì)分離
應(yīng)用假設(shè)條件:地震反射波和散射波具有相同的本質(zhì),但表現(xiàn)特征有異,它們是對(duì)同一入射波場(chǎng)在不同介質(zhì)尺度下地震波場(chǎng)響應(yīng)的不同表征[10,14-16]。在某種地震尺度上,地震反射波和散射波波場(chǎng)可相對(duì)分離。
地震散射波為地下非均質(zhì)體(異常擾動(dòng))的地震響應(yīng),它由多種類(lèi)型的地震波干涉疊加形成,信噪比較低,能量和頻率與反射波波場(chǎng)差異較大,通常散射波能量較弱,難以形成有效的強(qiáng)散射波場(chǎng),并淹沒(méi)于地震反射波剖面中,因此常常被忽略。而現(xiàn)行的地震波疊前成像方法具有類(lèi)似繞射波掃描疊加成像的功能,因此,其偏移成像結(jié)果已包含了豐富的散射波波場(chǎng)信息,它掩蓋于地震反射波成像剖面的背景能量中。擴(kuò)散濾波方法可通過(guò)調(diào)節(jié)擴(kuò)散系數(shù)和迭代次數(shù)實(shí)現(xiàn)地震反射波和散射波的相對(duì)分離。
圖2(a)為某地區(qū)經(jīng)地震保幅處理的疊前深度偏移剖面,圖2(b)為圖2(a)經(jīng)擴(kuò)散濾波方法處理后的地震剖面(反射波剖面),圖2(c)為圖2(a)經(jīng)擴(kuò)散濾波方法處理后的殘差地震剖面(散射波剖面)。對(duì)比擴(kuò)散濾波方法處理前和處理后的剖面(圖2)可知,擴(kuò)散濾波處理后的剖面更平滑[圖2(b)],但整體構(gòu)造特征和同相軸的變化方向與擴(kuò)散濾波前的剖面一致[圖2(a)],反映了該地區(qū)反射波波場(chǎng)的能量分布情況。與此同時(shí),反射波剖面的信噪比較處理前得到明顯改善,同相軸的連續(xù)性更好。而經(jīng)擴(kuò)散濾波處理后的殘差數(shù)據(jù)細(xì)節(jié)特征更清晰[圖2(c)],但其信噪比較反射波剖面更低,此雜亂反射的地震波能量主要為散射波能量,主要由該區(qū)的非均質(zhì)體引起。同時(shí),散射波剖面的垂向和橫向分辨率更高,斷層刻畫(huà)更清晰(圖2中的藍(lán)色圈),能夠有效地刻畫(huà)掩蓋于反射波能量中的地層散射點(diǎn)等精細(xì)地質(zhì)構(gòu)造信息(圖2中的紅色圈)。此外,地震散射波剖面上同相軸的抖動(dòng)位置正好代表地質(zhì)異常體和不規(guī)則地質(zhì)點(diǎn)的實(shí)際位置,符合地震反射界面不光滑的實(shí)際地質(zhì)與地球物理情況。
圖2 擴(kuò)散濾波方法在反射波和散射波波場(chǎng)分離中的應(yīng)用Fig.2 App lication of diffusion filteringmethod to wave field separation of reflection w ave and refraction wave
2.3 隨機(jī)噪音的有效壓制
應(yīng)用假設(shè)條件:隨機(jī)噪音是由地震采集或處理等過(guò)程引入到地震處理結(jié)果中的。
隨機(jī)噪音的存在直接降低了地震同相軸的連續(xù)性,影響后續(xù)的地震資料解釋等工作[17]。該噪音主要表現(xiàn)為豎條狀(采集過(guò)程引起,如采集腳印等),同相軸主要表現(xiàn)為鋸齒狀(主要由處理過(guò)程引起,如帶通濾波等)。擴(kuò)散濾波方法可通過(guò)調(diào)節(jié)擴(kuò)散系數(shù)和迭代次數(shù)實(shí)現(xiàn)保持?jǐn)鄬舆吘壍碾S機(jī)噪音壓制和同相軸連續(xù)性的有效增強(qiáng)。
圖3(a)為某地區(qū)的地震疊前時(shí)間偏移剖面,圖3(b)為擴(kuò)散濾波方法處理后的地震剖面,圖3(c)為擴(kuò)散濾波方法處理后的殘差地震剖面。分析圖3(a)可知,粉紅色橢圓內(nèi)存在豎條狀的噪音,同相軸橫向連續(xù)性較差,該種現(xiàn)象主要是由地震采集和數(shù)據(jù)處理過(guò)程中累積到最終成像結(jié)果的干擾能量引起的,常規(guī)處理方法較難將其壓制。采用各向異性擴(kuò)散濾波方法處理后,斷層兩側(cè)和連續(xù)同相軸兩側(cè)自動(dòng)采用兩種擴(kuò)散系數(shù),即沿邊緣方向采用較大的擴(kuò)散系數(shù),沿垂直邊緣方向采用較小的擴(kuò)散系數(shù),從而保證斷層邊緣得到有效保持[圖3(b)],連續(xù)的同相軸能量未被削弱,噪音得到有效壓制。且經(jīng)其處理的殘差剖面中均為地震隨機(jī)噪音[圖3(c)],無(wú)明顯的斷面和連續(xù)同相軸。由此可見(jiàn),各向異性擴(kuò)散濾波方法在地震噪音壓制中可取得更好的效果。
圖3 擴(kuò)散濾波方法在隨機(jī)噪音壓制中的應(yīng)用Fig.3 Application of diffusion filteringmethod to random noise suppression
2.4 地震數(shù)據(jù)插值
應(yīng)用假設(shè)條件:受地形起伏和障礙物等施工條件的綜合影響,在數(shù)據(jù)采集中常在某些位置出現(xiàn)地震記錄的道間距過(guò)大、異常道、數(shù)據(jù)不完整或數(shù)據(jù)缺失等問(wèn)題,這意味著一些地震響應(yīng)信息的丟失,并可能在各種處理過(guò)程中產(chǎn)生許多不必要的噪音。從目前地震插值技術(shù)的應(yīng)用現(xiàn)狀來(lái)看,地震道插值方法可在一定程度上恢復(fù)該地震波場(chǎng)信息,從而可提高偏移成像的精度,并切實(shí)有效地壓制噪音和抑制空間假頻出現(xiàn)[18]。而擴(kuò)散濾波方法在圖像處理領(lǐng)域得到了廣泛的應(yīng)用,其中一項(xiàng)為圖像修復(fù)技術(shù)[19-20]。筆者首次將該思路引入到地震處理領(lǐng)域,即對(duì)缺失的數(shù)據(jù)通過(guò)擴(kuò)散濾波方法處理實(shí)現(xiàn)修復(fù),從而在后續(xù)地震處理中避免假頻等問(wèn)題的出現(xiàn)。
該種基于圖像修復(fù)地震插值的基本原理可表示為:根據(jù)待修復(fù)區(qū)域周?chē)词軗p的地震數(shù)據(jù)來(lái)修正待修復(fù)區(qū)域中的每一個(gè)地震采樣點(diǎn)數(shù)值,且修復(fù)過(guò)程只改變待修復(fù)區(qū)域內(nèi)的地震數(shù)值,未受損區(qū)域內(nèi)的地震數(shù)值保持不變。基于擴(kuò)散濾波方法的地震數(shù)據(jù)修復(fù)技術(shù)與熱傳導(dǎo)方程中能量的傳遞過(guò)程類(lèi)似,地震數(shù)據(jù)為在一定約束條件下以振幅能量的形式向待修復(fù)區(qū)域進(jìn)行擴(kuò)散,最終使整個(gè)地震數(shù)據(jù)達(dá)到穩(wěn)定狀態(tài)。該種約束條件是根據(jù)地震數(shù)據(jù)本身的結(jié)構(gòu)特性(擴(kuò)散函數(shù))來(lái)控制各空間方向的擴(kuò)散強(qiáng)度,然后將最新迭代得到的地震數(shù)據(jù)和前一次迭代得到的地震數(shù)據(jù)之間的差值與給定的閾值進(jìn)行比較,如果該差值小于給定的閾值,則可結(jié)束迭代,并得到修復(fù)處理后的地震數(shù)據(jù)。若該差值大于給定的閾值,則繼續(xù)迭代計(jì)算,直至小于給定的閾值后結(jié)束迭代,最終實(shí)現(xiàn)地震數(shù)據(jù)的插值處理。
圖4(a)為某地區(qū)的地震疊前時(shí)間偏移剖面和理論合成的單炮模擬記錄,圖4(b)為圖4(a)經(jīng)數(shù)據(jù)抽稀后的結(jié)果,圖4(c)為圖4(b)經(jīng)擴(kuò)散濾波方法插值后的結(jié)果。對(duì)比分析圖4(a)和圖4(c)可知,原始數(shù)據(jù)和經(jīng)擴(kuò)散濾波方法插值后的地震同相軸連續(xù)性等特征基本保持一致,由此可見(jiàn),擴(kuò)散濾波方法在疊前和疊后地震數(shù)據(jù)插值處理中具有較好的應(yīng)用價(jià)值。
圖4 擴(kuò)散濾波方法在地震數(shù)據(jù)插值中的應(yīng)用Fig.4 App lication of diffusion filteringmethod to seism ic data interpolation
(1)擴(kuò)散濾波方法可有效實(shí)現(xiàn)低頻逆時(shí)偏移噪音的壓制處理,恢復(fù)掩蓋于低頻背景能量中的有效地層細(xì)節(jié)特征。
(2)擴(kuò)散濾波方法可實(shí)現(xiàn)地震反射波和散射波波場(chǎng)的相對(duì)分離處理,從而獲得掩蓋于反射波能量中的散射波能量,提高對(duì)地層細(xì)節(jié)的刻畫(huà)。
(3)擴(kuò)散濾波方法可實(shí)現(xiàn)保斷層邊緣的隨機(jī)噪音壓制處理,提高地震同相軸的連續(xù)性,壓制由采集或處理等過(guò)程引入的隨機(jī)干擾能量。
(4)擴(kuò)散濾波方法可實(shí)現(xiàn)地震數(shù)據(jù)的插值處理,并給出地震數(shù)據(jù)修復(fù)的基本原理,在一定程度上修復(fù)缺失的地震波場(chǎng)信息,實(shí)現(xiàn)地震數(shù)據(jù)的插值。
[1]王麗坤.三維各向異性擴(kuò)散濾波方法研究與應(yīng)用[D].成都:成都理工大學(xué),2011.
[2]Lu Yanhong,Lu Wenkai.Edge-preserving polynomial fittingmethod tosuppress random seismic noise[J].Societyof Exploration Geophysics,2009,74(4):69-73.
[3]Perona P,Malik J.Scale-space and edge detection using anisotropic diffusion[J].IEEE Transactions on Pattern Analysis and Machine Intelligence,1990,12(7):629-639.
[4]Fehmers GC,Hocker CFW.Faststructural interpretationwith structure-oriented filtering[J].Geophysics,2003,68(4):1286-1293.
[5]LavialleO,Pop S,Germain C,etal.Seismic faultpreservingdiffusion
:[J].JournalofApplied Geophysics,2007,61(2):132-141.
[6]Kadlec B J,Dorn GA,TufoHM.Confidenceand curvature-guided levelsets for channelsegmentation[G].SEGAnnual International Meeting Expanded Abstracts,2008:879-883.
[7]孫夕平,杜世通,湯磊.相干增強(qiáng)各向異性擴(kuò)散濾波技術(shù)[J].石油地球物理勘探,2004,39(6):651-655.
[8]王緒松,楊長(zhǎng)春.對(duì)地震圖像進(jìn)行保邊濾波的非線(xiàn)性各向異性擴(kuò)散算法[J].地球物理學(xué)進(jìn)展,2006,21(2):452-457.
[9]張爾華,王偉,高靜懷,等.非線(xiàn)性各向異性擴(kuò)散濾波器用于三維地震資料噪聲衰減與結(jié)構(gòu)特征增強(qiáng)[J].地球物理學(xué)進(jìn)展,2010,25(3):866-870.
[10]陳可洋,楊微,吳清嶺,等.地震反射波與散射波波場(chǎng)分離方法初探[J].巖性油氣藏,2013,25(2):76-81.
[11]張強(qiáng),張毅.一種用修正各向異性擴(kuò)散方程抑制SAR相干斑的算法[J].中國(guó)科學(xué)院研究生院學(xué)報(bào),2011,28(2):223-227.
[12]陳可洋,吳清嶺,范興才,等.地震波疊前逆時(shí)偏移脈沖響應(yīng)研究與應(yīng)用[J].石油物探,2013,52(2):163-170.
[13]陳可洋.基于拉普拉斯算子的疊前逆時(shí)噪聲壓制方法[J].巖性油氣藏,2011,23(5):87-95.
[14]尹軍杰.地震散射波場(chǎng)特征的數(shù)值模擬研究[D].北京:中國(guó)地質(zhì)大學(xué),2005.
[15]劉鐵華.地震散射波理論與數(shù)值模擬[D].西安:長(zhǎng)安大學(xué),2010.
[16]陳可洋.井間彈性波波場(chǎng)散射特征數(shù)值模擬分析[J].巖性油氣藏,2011,23(3):91-96.
[17]陳可洋,吳清嶺,李來(lái)林,等.松遼盆地三維地震資料連片處理關(guān)鍵技術(shù)及其應(yīng)用效果分析[J].巖性油氣藏,2012,24(2):87-91.
[18]陳可洋,范興才,吳清嶺,等.提高逆時(shí)偏移成像精度的疊前插值處理研究與應(yīng)用[J].石油物探,2013,52(4):409-416.
[19]陳可洋,范興才,吳清嶺,等.基于擴(kuò)散濾波的多尺度分解和重構(gòu)方法及其應(yīng)用初探[J].計(jì)算物理,2013,30(6):855-861.
[20]劉庚農(nóng).圖像修復(fù)技術(shù)研究[D].合肥:合肥工業(yè)大學(xué),2011.
(本文編輯:李在光)
App lication of diffusion filteringmethod to the seism ic data processing
CHEN Keyang1,WU Peixi2,YANGWei2
(1.Research Institute of Exploration and Development,PetroChina Daqing Oilfield Co.Ltd.,Daqing 163712,Heilongjiang,China;
2.No.6Oil Production Plant,PetroChina DaqingOilfield Co.Ltd.,Daqing163712,Heilongjiang,China)
Based on anisotropic diffusion filteringequation,thispaperpresented the correspondingnumericalequation, together with its simplest anisotropic diffusion filtering numerical equation.Then,we applied diffusion filtering method to four aspects of application,such asmulti-scale low-frequency reverse-timemigration noise suppression, seismic reflection wave and scattered wave relatively separating processing,holding fault boundary random noise suppression and pre-stack and post-stack seismic data interpolation processing,and also presented each assumed condition of the applications.The computation resultsshow that the diffusion filteringmethod can be used to obtain good processing resultsunder differentassumed conditionsand instruct thepracticalseismic datadigitalprocessing.
diffusion filtering;noisesuppression;seismicdatainterpolation;wave fieldseparation;seismicexploration
P631.4 < class="emphasis_bold">文獻(xiàn)標(biāo)志碼:A
A
1673-8926(2014)01-0117-06
2013-08-15;
2013-11-25
國(guó)家重點(diǎn)基礎(chǔ)研究發(fā)展計(jì)劃(973)項(xiàng)目“火山巖油氣藏的形成機(jī)制與分布規(guī)律”(編號(hào):2009CB219307)資助
陳可洋(1983-),男,碩士,工程師,主要從事高精度地震波傳播模擬與逆時(shí)成像、多線(xiàn)程并行計(jì)算與模塊開(kāi)發(fā)、實(shí)際地震資料數(shù)字處理方法研究與應(yīng)用等方面的工作。地址:(163712)黑龍江省大慶市讓胡路區(qū)大慶石油勘探開(kāi)發(fā)研究院地震處理二室。電話(huà):(0459)5508524。E-mail:keyangchen@163.com。