全波形反演中最優(yōu)化傳輸函數(shù)研究

胡光輝,賀偉光

(中石化石油物探技術(shù)研究院有限公司,江蘇南京211103)

全波形反演(FWI)將疊前地震數(shù)據(jù)信息映射為地下結(jié)構(gòu)或介質(zhì)參數(shù)分布,是目前地震學(xué)中最強(qiáng)大的建模成像工具之一。在勘探地球物理學(xué)中,早在1975年就建立了從地面地震數(shù)據(jù)推斷內(nèi)部結(jié)構(gòu)的概念[1-2],但FWI的理論建立者仍被認(rèn)為是LAILLY[3]和TARANTOLA[4]。對(duì)于全波形反演一詞的理解存在兩個(gè)觀點(diǎn):第1個(gè)觀點(diǎn)堅(jiān)持認(rèn)為只有利用整個(gè)地震炮集數(shù)據(jù)的波形信息才能稱為全波形反演(full指的是整個(gè)地震炮集數(shù)據(jù));第2個(gè)觀點(diǎn)更寬泛一些,把所選擇的地震事件信息完全利用即可成為全波形反演(full暗含的意思是fully)。無(wú)論是哪種觀點(diǎn),都認(rèn)為波的傳播需要經(jīng)過(guò)求解波動(dòng)方程來(lái)模擬,利用其它退化的模擬方法(比如射線方程法、單程波方程法)都不能稱為全波形反演。在全波形反演的早期發(fā)展中,由于計(jì)算能力有限以及缺乏低頻和長(zhǎng)偏移距觀測(cè)數(shù)據(jù),因而其發(fā)展十分緩慢。自2000年以來(lái),由于全波形反演理論上的發(fā)展和計(jì)算能力的進(jìn)步,該方法逐步發(fā)展成為現(xiàn)代工業(yè)生產(chǎn)中的常規(guī)地震資料處理方法之一,尤其是海上地震資料。然而,無(wú)論是在全波形反演的早期發(fā)展階段還是如今的工業(yè)應(yīng)用階段,周跳(cycle skipping)問(wèn)題依然沒(méi)有完全解決[5-6]。盡管有一些精巧的反演策略[7](頻率逐漸增加和偏移距逐漸增加),但是傳統(tǒng)的最小二乘函數(shù)(L2函數(shù))在許多情況下依然會(huì)陷入局部極值。

從數(shù)學(xué)上講,全波形反演是在彈性動(dòng)力學(xué)方程約束下的局部?jī)?yōu)化問(wèn)題。數(shù)學(xué)上的處理有兩種:第1種思路是構(gòu)建拉格朗日函數(shù);第2種思路是建立增廣拉格朗日函數(shù)。增廣拉格朗日函數(shù)是將彈性動(dòng)力學(xué)方程作為二次項(xiàng)來(lái)構(gòu)建優(yōu)化函數(shù)[8]。從數(shù)學(xué)界來(lái)看,增廣拉格朗日函數(shù)是凸性良好的優(yōu)化函數(shù)。波動(dòng)方程以二次型的方式引入增廣拉格朗日函數(shù)并不要求波動(dòng)方程被精確地滿足。這種處理方式相當(dāng)于增大了波形反演的解的搜索空間,同時(shí)增強(qiáng)了方程的凸性。但是,增廣拉格朗日函數(shù)的計(jì)算要求很高,原因是涉及模擬算子和伴隨算子的乘積。在頻率域,計(jì)算代價(jià)可以接受,但在時(shí)間域則相當(dāng)困難。因此,增廣拉格朗日函數(shù)在工業(yè)生產(chǎn)中很少使用。全波形反演的主流框架一直是第1個(gè)思路,PLESSIX[9]對(duì)該框架進(jìn)行了綜述。在拉格朗日框架中,作為約束的彈性動(dòng)力學(xué)方程被線性地引入到優(yōu)化函數(shù)中,之后可以方便地得到伴隨場(chǎng)震源和梯度的計(jì)算。拉格朗日框架要求在每次迭代時(shí)都能準(zhǔn)確地求解彈性動(dòng)力學(xué)方程。等價(jià)地說(shuō),在全波形反演的整個(gè)迭代過(guò)程中,約束方程(彈性動(dòng)力學(xué)方程)始終得到滿足。這限制了未知參數(shù)的搜索空間,因此導(dǎo)致拉格朗日框架的凸性不如增廣拉格朗日框架。但是,拉格朗日框架因?yàn)橛?jì)算代價(jià)較低而成為學(xué)術(shù)界和工業(yè)界的主流方式。在拉格朗日函數(shù)框架方案中,避免周跳的主要方法是設(shè)計(jì)凸目標(biāo)函數(shù),如L1范數(shù)[10]、混合或Huber函數(shù)[11-13]、零延遲互相關(guān)[14-17](相當(dāng)于歸一化的L2函數(shù))、包絡(luò)函數(shù)[18-24]、懲罰互相關(guān)函數(shù)[25]、基于反褶積的函數(shù)[26]和最優(yōu)傳輸函數(shù)[27-33]。其中,最優(yōu)化傳輸函數(shù)的優(yōu)勢(shì)是其比較的是全局地震數(shù)據(jù),而不是逐點(diǎn)比較地震數(shù)據(jù)。這種特性保持了地震數(shù)據(jù)的事件間和檢波點(diǎn)端記錄數(shù)據(jù)的相干性,從而使得較傳統(tǒng)L2函數(shù)能更穩(wěn)定地反演。本研究采用拉格朗日函數(shù)的框架,綜述了幾種最優(yōu)傳輸函數(shù)在全波形反演中的應(yīng)用。

最優(yōu)化傳輸起源于法國(guó)的一個(gè)工程問(wèn)題,即如何將沙子從源頭位置運(yùn)輸出來(lái),以最少的能耗來(lái)填補(bǔ)施工區(qū)域的洞。這個(gè)問(wèn)題被MONGE[34]用數(shù)學(xué)方法表述為一個(gè)L1最小化問(wèn)題。在一些研究中,它可能被命名為EMD(Earth Mover Distance)。方程的答案可以通過(guò)求解非線性方程,這種非線性方程也稱為Monge-Ampere方程。將1范數(shù)推廣到任何其它范數(shù)稱為p-Wasserstein函數(shù),其中,p表示范數(shù)。證明Monge方程解的存在性,以及求解都是很困難的。KANTOROVICH[35]取得了突破,將Monge公式擴(kuò)展到更高維空間會(huì)得到線性問(wèn)題;具體來(lái)說(shuō),這是線性分配問(wèn)題。FERRADANS等[36]和LELLMANN等[37]在圖像處理中引入最優(yōu)化傳輸函數(shù),該函數(shù)是一種自動(dòng)且正確的匹配模式。最優(yōu)化傳輸函數(shù)目前已廣泛應(yīng)用于機(jī)器學(xué)習(xí)或人工智能(AI)等領(lǐng)域。最優(yōu)傳輸問(wèn)題根據(jù)分布是連續(xù)的還是離散的[38]分為3類,即連續(xù)(兩個(gè)分布都是連續(xù)的)、半離散(一個(gè)是連續(xù)的,另一個(gè)是離散的)和離散的(兩者都是離散的)。連續(xù)問(wèn)題可以用Monge-Ampere方程或流體動(dòng)力學(xué)方程來(lái)解決。半離散問(wèn)題利用計(jì)算幾何和牛頓優(yōu)化技術(shù)可以快速解決。離散問(wèn)題用正則化最優(yōu)傳輸求解器或標(biāo)準(zhǔn)線性賦值求解器求解。

直接計(jì)算最優(yōu)傳輸問(wèn)題的運(yùn)算量巨大。例如,當(dāng)將n個(gè)點(diǎn)移動(dòng)到n個(gè)目標(biāo)位置時(shí),每個(gè)點(diǎn)都可以移動(dòng)到任何位置。因此,未知變量的數(shù)目為n2。KANTOROVICH[35]將原始的最優(yōu)化傳輸問(wèn)題轉(zhuǎn)換到對(duì)偶空間中,未知量的個(gè)數(shù)變?yōu)閚,極大地減少了離散情形下未知變量的個(gè)數(shù)。在對(duì)偶空間中,這可以用標(biāo)準(zhǔn)的線性方程求解器來(lái)方便地求解。同時(shí),也有一些專門針對(duì)線性分配的算法加速收斂,比如匈牙利算法[39]或拍賣算法[40]。在圖像比較領(lǐng)域,在傳輸代價(jià)函數(shù)中添加關(guān)于傳輸算子的正則化被廣泛采用[41]。為了提高計(jì)算效率,人們發(fā)展了更為精巧和復(fù)雜的計(jì)算技術(shù),如多層次采樣策略、多尺度正則化策略。在多層次采樣策略中,將整個(gè)數(shù)據(jù)集采樣到多個(gè)層次。將原始數(shù)據(jù)集定義為零級(jí),沿每個(gè)維度以2的比率進(jìn)行降采樣得到第1級(jí)(對(duì)于2D問(wèn)題,這個(gè)數(shù)字現(xiàn)在是原始數(shù)據(jù)集的0.25倍)。第1級(jí)再次重新采樣以獲得下一級(jí),最后一級(jí)是最小的數(shù)據(jù)集合。數(shù)據(jù)分級(jí)采樣之后,求解器從最小的數(shù)據(jù)集開(kāi)始(也是最終的數(shù)據(jù)集)。最后一級(jí)的傳輸平面計(jì)算出來(lái)之后,進(jìn)行插值,成為前一級(jí)的初始會(huì)傳輸平面。這個(gè)過(guò)程一直持續(xù)到零級(jí)。在每個(gè)層面上,都采用多尺度的正則化策略來(lái)計(jì)算傳輸平面。在各級(jí)計(jì)算中,正則化方法相應(yīng)地從相對(duì)較大的值開(kāi)始,逐漸減小到接近于零。在傳輸平面代表的映射中,每個(gè)點(diǎn)會(huì)被分割到幾個(gè)位置上。在零級(jí)數(shù)據(jù)集上,傳輸平面上會(huì)有一些非常集中的能量分布,每個(gè)能量團(tuán)的中心對(duì)應(yīng)于一個(gè)一一映射的解。正則化的引入使得無(wú)法得到單點(diǎn)對(duì)單點(diǎn)的映射,但是可以定義為質(zhì)心坐標(biāo)強(qiáng)制傳輸平面成為一對(duì)一的映射。

最優(yōu)化傳輸并不能直接應(yīng)用到地震數(shù)據(jù)處理中。地震信號(hào)是上、下震蕩,有正、有負(fù)的。而最優(yōu)化傳輸處理的是正信號(hào),且兩個(gè)比較對(duì)象的質(zhì)量要相等。因此需要對(duì)數(shù)據(jù)進(jìn)行預(yù)處理,使得地震信號(hào)滿足非負(fù)且質(zhì)量相等的條件。依據(jù)數(shù)據(jù)預(yù)處理的方式,目前有3種最優(yōu)化傳輸公式在全波形反演中得到了很好的應(yīng)用,即2-Wasserstein函數(shù)[27-29],KR-OT函數(shù)[30-32]和圖空間的OT函數(shù)GSOT[33]。2-Wasserstein函數(shù)是基于范數(shù)2準(zhǔn)則的最優(yōu)傳輸函數(shù)的公式。它是第一個(gè)被引入到地球物理學(xué)的最優(yōu)傳輸函數(shù)[27](W22函數(shù),第一個(gè)2意味著2范數(shù),第二個(gè)2意味著平方)。它是一維的算法,每次對(duì)比一對(duì)地震道信號(hào)。每一個(gè)合成地震道和觀測(cè)地震道都應(yīng)該經(jīng)過(guò)預(yù)處理為正的并歸一化(保證質(zhì)量相等)。最簡(jiǎn)單的預(yù)處理是減去最小值并除以其積分。該公式的獨(dú)特優(yōu)點(diǎn)是存在解析解,計(jì)算速度快且穩(wěn)定。擴(kuò)展到其它規(guī)范(p-Wasserstein函數(shù))是很直接的。缺點(diǎn)是預(yù)處理對(duì)理論凸性略有破壞。在極端情形下,可以假設(shè)地震信號(hào)減去負(fù)無(wú)窮大,在此情形下,最優(yōu)化傳輸函數(shù)退化為常規(guī)的最小二乘誤差。

第2個(gè)引入地球物理學(xué)的最優(yōu)化傳輸是由METIVIER等[30]提出的KR-OT。它是基于1范數(shù)(1-Wasserstein函數(shù))的最優(yōu)化傳輸函數(shù),并且是在對(duì)偶空間中表示。KR-OT在一維、二維和三維空間中都具有相似的數(shù)學(xué)結(jié)構(gòu)。二維全波形反演的案例表明,2D KR-OT明顯優(yōu)于1D KR-OT。正如METIVIER所解釋的那樣,1D KR-OT反演的退化是忽略了接收器間的相干性。在KR-OT的公式中,最優(yōu)化傳輸函數(shù)只依賴于兩個(gè)數(shù)據(jù)殘差(地震數(shù)據(jù)可以直接發(fā)送到求解器,而不需要移動(dòng)或歸一化),這是優(yōu)點(diǎn)也是缺點(diǎn)。弊端在于所有的地震事件從一開(kāi)始就混合在一起,當(dāng)初始合成數(shù)據(jù)與觀測(cè)數(shù)據(jù)相差較遠(yuǎn)時(shí),映射可能無(wú)法正確匹配地震事件。KR-OT的局限之一是采樣在時(shí)間軸和空間軸上應(yīng)該是均勻的。至于將KR-OT擴(kuò)展到其它領(lǐng)域,如超聲成像[42],原則上是不合適的,因?yàn)榻邮掌鞑⒉皇蔷鶆蚍植嫉摹?/p>

上述兩個(gè)最優(yōu)化傳輸公式都屬于連續(xù)傳輸問(wèn)題。勘探地球物理學(xué)中的第3種最優(yōu)傳輸函數(shù)是圖空間最優(yōu)傳輸函數(shù)GS-OT,這是離散傳輸問(wèn)題。將一維離散信號(hào)擴(kuò)展為時(shí)間-振幅二維圖中的點(diǎn)集[33,43-44]。這種處理方法自然滿足了數(shù)據(jù)非負(fù)的要求。同時(shí),如果用相同數(shù)量的點(diǎn)離散地震數(shù)據(jù),質(zhì)量守恒也是自然滿足的。從原理上將沿著任何一個(gè)維度都可以任意采樣,去除了均勻采樣的要求。GS-OT僅取決于點(diǎn)的數(shù)量,即在一個(gè)地震道中的離散采樣點(diǎn)。所以直接擴(kuò)展到高維空間必然面對(duì)計(jì)算量巨大的挑戰(zhàn)。

最近,XU等[45]開(kāi)發(fā)了基于圖空間公式的高維最優(yōu)傳輸函數(shù)。雖然這項(xiàng)研究是在全波形反演中的應(yīng)用,但作者最初的研究動(dòng)機(jī)是要匹配兩個(gè)彩色點(diǎn)云。每個(gè)顏色點(diǎn)都與一個(gè)向量(R,G,B)相關(guān)聯(lián),但最優(yōu)化傳輸只針對(duì)標(biāo)量函數(shù)。從數(shù)學(xué)上講,XU等[45]開(kāi)發(fā)的是比較兩個(gè)向量場(chǎng)的最優(yōu)化傳輸函數(shù)。該研究進(jìn)行了二維和三維實(shí)驗(yàn)。正如作者解釋的那樣,當(dāng)點(diǎn)的數(shù)量少于1000時(shí),計(jì)算非?？?但對(duì)于超過(guò)10000個(gè)點(diǎn)的彩色點(diǎn)云來(lái)說(shuō),計(jì)算時(shí)間過(guò)長(zhǎng)。因?yàn)樵撗芯繉⒄麄€(gè)地震圖分成互相不重疊的小區(qū)域(二維地震圖是長(zhǎng)方形,三維是長(zhǎng)方體),大幅度降低了計(jì)算成本。依據(jù)在于,地震圖的對(duì)比中,合成數(shù)據(jù)的地震事件是不應(yīng)當(dāng)與觀測(cè)事件相差太遠(yuǎn)。所以劃成小區(qū)域后,數(shù)據(jù)匹配依然存在合理性。但是在圖像對(duì)比中,整張圖可能會(huì)旋轉(zhuǎn)或者倒轉(zhuǎn),在這種情形下不能劃分小區(qū)域。經(jīng)過(guò)數(shù)據(jù)區(qū)域劃分處理后,將高維最優(yōu)傳輸函數(shù)納入全波形反演中變得可行和現(xiàn)實(shí)。該研究還將熵變換應(yīng)用于空間坐標(biāo)和時(shí)間坐標(biāo)。雖然在一些模型中是否采用熵變換并不影響反演,但有趣的是,在Chevron 2014盲測(cè)中,應(yīng)用熵變換的反演結(jié)果優(yōu)于沒(méi)有熵變換的反演。

最優(yōu)化傳輸函數(shù)是在原始空間或?qū)ε伎臻g中求解。在數(shù)學(xué)領(lǐng)域,一種新的基于混合算法的最優(yōu)傳輸方法被提出。該方法是基于2-Wasserstein函數(shù)的二維算法,不再是一維單道處理。這是由JACOBS等[46]開(kāi)發(fā)的一種新穎的方法來(lái)解決最優(yōu)化傳輸函數(shù)。首先,將原始形式轉(zhuǎn)換到對(duì)偶空間,在數(shù)學(xué)上就得到了兩個(gè)等價(jià)的方程:一個(gè)在原始空間,另一個(gè)在對(duì)偶空間。新穎的部分在于第2步。JACOBS等[46]只迭代一步原始解;然后跳到對(duì)偶空間迭代一步;之后再回到原始空間繼續(xù)循環(huán),直到得到一個(gè)收斂的結(jié)果。JACOBS等[46]將其命名為往復(fù)法(back and forth method,BFM)。

最初的最小二乘目標(biāo)函數(shù)是逐點(diǎn)比較的算法,后來(lái)興起的從互相關(guān)到最優(yōu)傳輸?shù)哪繕?biāo)函數(shù)更多是基于全局比較算法。因此,可以更好地利用相鄰點(diǎn)之間甚至地震事件之間的相干性。這很好地提高了目標(biāo)函數(shù)的凸性和地下結(jié)構(gòu)的一致性。本研究旨在提供對(duì)全波形反演中最優(yōu)化傳輸函數(shù)的介紹,解釋了目前引入全波形反演的幾種最優(yōu)化傳輸公式,包括W2(1D版本,2-Wasserstein函數(shù))、BFM、KR-OT、GS-1D和GS-2D(高維最優(yōu)傳輸函數(shù))。

1 方法理論

本文重點(diǎn)在于討論最優(yōu)化傳輸函數(shù),但為了說(shuō)明與全波形反演的結(jié)合應(yīng)用,先簡(jiǎn)單介紹全波形反演的實(shí)現(xiàn)。之后重點(diǎn)介紹5個(gè)最優(yōu)化傳輸函數(shù)的基本原理。W2是1D 2-Wasserstein函數(shù),BFM是2D或者3D的2-Wasserstein函數(shù),KR-OT是2D或者3D的1-Wasserstein函數(shù),GS-1D是一維信號(hào)的圖空間最優(yōu)化傳輸,GS-2D是多維信號(hào)的圖空間最優(yōu)化傳輸。

1.1 全波形反演

FWI已經(jīng)發(fā)展了幾十年,現(xiàn)在是地震資料處理業(yè)界的常規(guī)技術(shù)。FWI包含一個(gè)目標(biāo)函數(shù)和一個(gè)約束方程(又稱為控制方程)。目標(biāo)函數(shù)是刻畫合成地震數(shù)據(jù)與觀測(cè)地震數(shù)據(jù)直接偏差的函數(shù),控制方程是描述聲波或彈性地震波傳播的方程。最早也是最常見(jiàn)的目標(biāo)函數(shù)是最小二乘函數(shù)(L2函數(shù)),即:

(1)

式中:p(t)是合成數(shù)據(jù);d(t)是觀測(cè)數(shù)據(jù),求和針對(duì)所有的震源和接收器。目標(biāo)函數(shù)中應(yīng)用了加權(quán)函數(shù)W(t),因?yàn)椴⒉皇撬械牡卣鹗录急焕?。例?在應(yīng)用聲波方程求解器處理陸地地震數(shù)據(jù)時(shí),需要預(yù)先去除觀測(cè)數(shù)據(jù)中的面波和橫波。目標(biāo)函數(shù)的優(yōu)化通常是基于局部?jī)?yōu)化求解器,模型更新的搜索方向是基于梯度的。梯度的計(jì)算通常是采用伴隨方[9]。正傳波場(chǎng)與伴隨波場(chǎng)之間的累積求和為梯度。計(jì)算梯度只需要花費(fèi)2～3個(gè)波傳播的計(jì)算時(shí)間。伴隨波場(chǎng)是通過(guò)反向傳播伴隨源,伴隨源由目標(biāo)函數(shù)關(guān)于合成數(shù)據(jù)的一階導(dǎo)數(shù)給出。目標(biāo)函數(shù)的一階擾動(dòng)為:

(2)

由此可確定伴隨源s(t)=W(t)[p(t)-d(t)]。確定了伴隨源之后,即可通過(guò)正傳波場(chǎng)與反傳波場(chǎng)的相關(guān)求取模型更新的梯度,繼而通過(guò)不同的迭代求解算法逐步更新模型。

1.2 最優(yōu)化傳輸函數(shù)的模式匹配特性

最優(yōu)傳輸函數(shù)基于模式的匹配。這個(gè)性質(zhì)可以用插值的例子來(lái)說(shuō)明。圖1在左下角顯示圓盤,在右上角顯示正方形盤。傳統(tǒng)的插值(線性插值)將它們與兩個(gè)加權(quán)因子分別相乘,之后相加。圖2中的第1行顯示了線性插值的演化,其中圓盤上的加權(quán)因子從1逐漸減小到0,同時(shí)正方形盤上的加權(quán)因子從0增加到1?？梢钥吹?插值的效果是此消彼長(zhǎng)的變化。但有時(shí),希望可以像水流一樣把物體形態(tài)光滑地變化成為另一個(gè)。圖2中的第2行顯示的正是這種變化,這一種基于最優(yōu)化傳輸?shù)牟逯?首先求出圓盤與正方形盤直接的傳輸平面,然后對(duì)平面插值即可?？梢钥吹?圓盤的形狀和位置均平滑地演化成正方形盤和其所在的位置。

最優(yōu)傳輸函數(shù)的這種插值特性對(duì)FWI非常有幫助。傳統(tǒng)的FWI經(jīng)常被周跳問(wèn)題困擾。最小二乘函數(shù)只是點(diǎn)對(duì)點(diǎn)地比較兩個(gè)地震道。一旦兩個(gè)地震事件相距較遠(yuǎn),在伴隨場(chǎng)震源上兩個(gè)事件就不重疊,也就是說(shuō),它們之間就沒(méi)有相互作用了。這就解釋了為什么最小二乘函數(shù)只在兩個(gè)地震事件誤差在半周期內(nèi)時(shí)起作用,同時(shí)在誤差大于半周期的時(shí)候會(huì)陷入局部極值。而波動(dòng)方程層析建模的凸性較FWI更好,因?yàn)閮蓚€(gè)事件總是可以通過(guò)互相關(guān)聯(lián)系起來(lái)。因此,我們期望具有全局對(duì)比功能的最優(yōu)化傳輸函數(shù),它比最小二乘函數(shù)具有更好的凸性。

1.2.1W2函數(shù)

(3)

應(yīng)遵守質(zhì)量守恒原則,即:

(4)

圖3 計(jì)算W2的示意

基于方程(3)和方程(4),可以計(jì)算出伴隨源。伴隨源s(t)由W2相對(duì)于合成的p(t)的偏導(dǎo)數(shù)給出[47]。引入中間變量,有:

(5)

最后的伴隨源為:

(6)

數(shù)值實(shí)現(xiàn)過(guò)程可總結(jié)為:

1) 對(duì)每個(gè)地震道進(jìn)行預(yù)處理,使其為正。建議將每個(gè)處理過(guò)的地震道除以其自身的質(zhì)量,使總質(zhì)量為1。

移動(dòng)和歸一化屬于線性變換,把數(shù)據(jù)變?yōu)檎盘?hào)還包括非線性變換。非線性變換包括平方、指數(shù)或雙曲正切(tanh)。CHEN等[47]在地震定位的研究中采用平方法對(duì)地震數(shù)據(jù)進(jìn)行預(yù)處理。作者設(shè)計(jì)了一個(gè)兩層模型,并根據(jù)模型層位速度的不同計(jì)算了相對(duì)應(yīng)的合成地震圖,最終做了目標(biāo)函數(shù)的凸性分析,如圖4 所示。令人驚訝的是,將地震數(shù)據(jù)進(jìn)行平方這個(gè)預(yù)處理優(yōu)于移動(dòng)和歸一化策略。采用平方預(yù)處理得到的目標(biāo)函數(shù)分布圖更加光滑,凸性更好。

圖4 目標(biāo)函數(shù)誤差分布[47]

與此同時(shí),體波的振幅和弱事件的振幅會(huì)得到增強(qiáng)。振幅平衡是彈性全波形反演中重要的步驟之一,對(duì)于提升體波對(duì)于深部的反演貢獻(xiàn)巨大。HE等[31-32]在各向異性介質(zhì)中詳細(xì)論述了振幅問(wèn)題。

1.2.2BFM函數(shù)

一維問(wèn)題的最優(yōu)化傳輸存在解析解。與一維問(wèn)題不同,高維問(wèn)題不存在最優(yōu)化傳輸?shù)慕馕鼋?。多年?lái),數(shù)學(xué)家們一直在努力尋找有效的求解方法。來(lái)回跳躍法BFM屬于有效求解的方法之一。首先,對(duì)原問(wèn)題構(gòu)造拉格朗日函數(shù),將其轉(zhuǎn)化為對(duì)偶空間的等價(jià)問(wèn)題。至此,有兩種等價(jià)的最優(yōu)化傳輸問(wèn)題,一種在原始空間,另一種在對(duì)偶空間。BFM的獨(dú)特優(yōu)點(diǎn)是,最優(yōu)解是在原始空間和對(duì)偶空間之間迭代產(chǎn)生的。在原始空間迭代一次,然后跳到對(duì)偶空間,在對(duì)偶空間同樣只迭代一次,然后再跳回原始空間;一直循環(huán)直到問(wèn)題收斂。這種方法的數(shù)學(xué)基礎(chǔ)比較復(fù)雜[46]。BFM仍然要求輸入是非負(fù)的,用于W2的預(yù)處理也適用于BFM。

1.2.3KR-OT函數(shù)

KR-OT是基于1-Wasserstein的函數(shù);但并不是求解原始問(wèn)題,而是轉(zhuǎn)化到對(duì)偶空間中求解。KR-OT的代價(jià)函數(shù)是基于1范數(shù)c(x,y)=|x-y|,在對(duì)偶空間中的表達(dá)式為[30]:

(7)

式中:φ(x)在有界的1-Lipschitz空間中,而λ是一個(gè)用戶定義的參數(shù)。第1個(gè)不等式保證了φ不能太大的要求,第2個(gè)不等式排除了突然變化的情況。HE等[32]仔細(xì)研究了不等式約束的影響,其結(jié)論是存在較大的λ取值空間,都可以給出合理的解決方案φ。KR-OT的顯著特性是,該公式只依賴于兩個(gè)分布之間的差異,所以不要求輸入的地震信號(hào)是正的。對(duì)于任何負(fù)分布或非正分布,可以將兩個(gè)分布加上同一個(gè)常數(shù)變?yōu)檎?但是從KR-OT的表達(dá)式看,這個(gè)常數(shù)是會(huì)被消去而不起任何作用。

在全波形反演中,地震數(shù)據(jù)是以離散形式存儲(chǔ)的。在三維情況下,對(duì)應(yīng)的KR-OT公式為:

(8)

式中:常數(shù)因子Δx1,Δx2,Δx3在公式(8)中的第1個(gè)方程中被省略,uijk是合成數(shù)據(jù)與觀測(cè)數(shù)據(jù)的差異。公式(8)可以用SDMM(simultaneous direction method of multipliers)高效求解。METIVIER等[30]詳細(xì)說(shuō)明了如何解決KR-OT。其中,求解矩陣逆是一個(gè)沉重的計(jì)算負(fù)擔(dān)。后來(lái),METIVIER等[48]證明了其中的矩陣運(yùn)算與泊松方程是等價(jià)的,從而可以在頻率域用FFT方法或在時(shí)域用多重網(wǎng)格方法有效地求解泊松方程。雖然如此,三維計(jì)算仍然非常昂貴,主要原因是每次求解KR-OT都需要多次求解三維泊松方程。在二維空間中,從公式(8)中去掉一維,可以大大減少計(jì)算時(shí)間。因此,偽三維處理也是折衷的方法。將三維地震數(shù)據(jù)體沿地震數(shù)據(jù)線方向分割成切片(地震數(shù)據(jù)切片形成二維矩陣),并使用2D KR-OT求解器對(duì)每個(gè)切片進(jìn)行處理。

KR-OT的輸入信號(hào)是殘差,故所有的地震事件從一開(kāi)始就是混合的。這會(huì)影響KR-OT的凸性。與其它OT公式相比,KR-OT的獨(dú)特優(yōu)勢(shì)是KR-OT將解約束為1-Lipschitz函數(shù)。這種優(yōu)勢(shì)對(duì)于包含面波的陸地地震數(shù)據(jù)具有很大的潛力。圖5為KR-OT能平衡振幅,展示的是一個(gè)陸地?cái)?shù)據(jù)全波形反演的伴隨場(chǎng)。L2伴隨源以面波為主,因此,重建由體波采樣的深層區(qū)域是困難的。相比之下,KR-OT伴隨源的振幅分布平衡,反演確實(shí)顯示了從淺到深的重構(gòu)模型。

圖5 KR-OT能平衡振幅[32]

1.2.4GS-1D函數(shù)

圖空間最優(yōu)化傳輸(Graph Space Optimal Transport)處理的是離散化的問(wèn)題。將一維地震道擴(kuò)展為時(shí)間-振幅的二維空間,最優(yōu)化傳輸?shù)拇鷥r(jià)函數(shù)為:

cij=(pi-dj)2+η(ti-tj)2

(9)

式中:pi是在ti時(shí)刻采樣的合成地震道;dj是在tj時(shí)刻采樣的觀測(cè)地震道;cij是點(diǎn)(ti,pi)至點(diǎn)(tj,dj)的移動(dòng)成本;η是用戶定義的參數(shù),用來(lái)控制時(shí)間軸的權(quán)重。最優(yōu)化傳輸要求質(zhì)量守恒,在圖空間中,質(zhì)量守恒約束指的是合成地震道與觀測(cè)地震道包含相同數(shù)量的點(diǎn)。配對(duì)一對(duì)地震道的總成本是對(duì)所有個(gè)體成本的求和,即:

(10)

(11)

伴隨源si看起來(lái)像傳統(tǒng)的L2的情況,但其中的觀測(cè)數(shù)據(jù)是重排列之后的觀測(cè)數(shù)據(jù)。在L2的情況下,兩個(gè)不重疊地震事件是沒(méi)有相互作用的,這也是周跳的根源。GS 1D依然能夠通過(guò)求解最優(yōu)的排列來(lái)聯(lián)系這些事件。

圖6顯示了L2、KR-OT和GS-1D的伴隨源。圖6a中,觀測(cè)數(shù)據(jù)包含多個(gè)地震事件,而合成數(shù)據(jù)只有兩個(gè)。圖6b中,常規(guī)的最小二乘函數(shù)無(wú)法平衡振幅,而KR-OT和GS-OT都在很好的平衡振幅。輸入信號(hào)由Ricker子波組成,模擬的是事件錯(cuò)位和事件缺失的情況。L2伴隨源具有較大的振幅差異,而OT伴隨源確實(shí)平衡了振幅。然而,這并不意味著在這兩個(gè)OT函數(shù)中都忽略了振幅信息。在全波形反演中,隨著迭代的繼續(xù),觀測(cè)數(shù)據(jù)的振幅會(huì)得到越來(lái)越好的利用。最后,所有的事件都是完全擬合的。

圖6 L2、KR-OT和GS-OT的伴隨源[31]

1.2.5GS-2D(高維OT)

(12)

圖7和圖8均引自XU 等[45]。圖7顯示了兩個(gè)二維質(zhì)量分布之間的映射;兩個(gè)質(zhì)量分布都是由高斯函數(shù)生成,但是中心位置不同。二維的GS-OT算法很好地求解出了匹配關(guān)系(黑色實(shí)線)。若是采用一維算法,則只能沿著橫向或者縱向匹配,而無(wú)法正確求解出真實(shí)的匹配關(guān)系。在全波形反演中,地震事件經(jīng)常是在時(shí)間軸方向和空間軸方向都有移動(dòng),這也是為何要發(fā)展高維OT算法的原因。圖8顯示的是三維質(zhì)量分布之間的三維映射。兩個(gè)球體分別對(duì)應(yīng)的是兩個(gè)三維的質(zhì)量分布,可以看出,三維的GS-OT算法正確地匹配到了兩個(gè)三維質(zhì)量分布。XU 等[45]將高維最優(yōu)化傳輸GS 2D應(yīng)用于聲學(xué)和彈性全波形反演試驗(yàn),包括彈性SEG/EAGE逆沖模型和Chev-ron 2014的數(shù)據(jù)盲測(cè)試驗(yàn)。在逆沖模型中,L2反演重構(gòu)出了橫波速度參數(shù)vS的淺層結(jié)構(gòu),GS 1D重建了vS的淺部和深部。相應(yīng)地,體波和面波都比L2得到了更好的利用。GS-2D得到了最佳的反演,即vP和vS都顯示了高分辨率和高保真度的結(jié)構(gòu)。在XU等[45]展示的Chevron 2014數(shù)據(jù)盲測(cè)實(shí)驗(yàn)中,使用熵變換的代價(jià)函數(shù)x=xlogx得到了更干凈的速度。

圖7 二維質(zhì)量分布之間的最優(yōu)化映射

圖8 三維質(zhì)量分布之間的最優(yōu)化映射

2 反演實(shí)驗(yàn)

我們已經(jīng)實(shí)現(xiàn)了上述所有的最優(yōu)化傳輸函數(shù)W2,BFM,KR-OT,GS-1D和GS-2D的FWI。本節(jié)將展示兩個(gè)模型測(cè)試結(jié)果。第1個(gè)是透射例子;第2個(gè)是選取的SEAM Ⅱ Foothill模型的切片。

2.1 透射模型

真實(shí)模型是在均勻背景vP=3.0km/s的模型中添加了vP=3.6km/s的圓形異常體。圓形異常體的直徑為3.6km。假設(shè)密度為常數(shù)ρ=1g/cm3。以均勻背景速度作為初始模型。圖9顯示了真實(shí)的模型和初始的模型。在頂部等距離布置了32個(gè)爆炸源,在底部布置了300個(gè)壓力接收器。震源子波為10Hz的Ricker子波。以10Hz為參考計(jì)算,圓形異常體的大小是主波長(zhǎng)的10倍。

圖9 全波形反演測(cè)試

圖10顯示了觀測(cè)數(shù)據(jù)以及觀測(cè)數(shù)據(jù)與合成數(shù)據(jù)的交錯(cuò)比較。對(duì)比圖清楚地說(shuō)明了周期跳躍問(wèn)題。L2反演進(jìn)行了6次迭代陷入了局部極值,反演結(jié)果在圖11a中。KR-OT反演也被局限在某個(gè)局部最小值中(圖11d)。在這個(gè)特定的例子中,初始模型與實(shí)際模型相差太遠(yuǎn)。正如作者所解釋的,KR-OT將L2誤差作為輸入,它的性能是受到影響的。其它OT反演(W2,BFM,GS-1D,GS-2D)都正確地恢復(fù)了模型。

圖10 觀測(cè)數(shù)據(jù)、合成數(shù)據(jù)與觀測(cè)數(shù)據(jù)對(duì)比

2.2 SEAM Ⅱ Foothill模型

所用模型是從原始的3D SEAM Ⅱ FOOTHILL模型中抽取的一個(gè)剖面,同時(shí)去除了近地表低速區(qū)。真實(shí)模型如圖12所示,縱、橫向采樣均是20m,橫向有720個(gè)格點(diǎn),長(zhǎng)約14km;縱向有200個(gè)格點(diǎn),深約4km?？梢悦黠@看到地層發(fā)生扭曲和位移,但是模型沿著縱向的整體變化率并不大。采用3Hz的Ricker子波作為震源子波。在頂部布置了48個(gè)震源和270個(gè)接收器。模型正演采用了吸收邊界條件。圖13為觀測(cè)數(shù)據(jù)的炮集。圖中的初至大致沿著一條直線,因?yàn)槟Ｐ涂v向變化不大,所以潛水波(diving wave)穿透深度較淺,從而導(dǎo)致走時(shí)跟距離接近線性關(guān)系。本文創(chuàng)建了4個(gè)初始模型,如圖14 所示。前3個(gè)由真實(shí)模型采用高斯平滑得到,平滑半徑分別為10,20,30個(gè)網(wǎng)格長(zhǎng)度。最后一個(gè)是一維線性遞增模型。相應(yīng)的合成地震數(shù)據(jù)如圖15所示。初始模型越簡(jiǎn)單,合成數(shù)據(jù)也越簡(jiǎn)單。對(duì)前面兩個(gè)初始模型,尚能識(shí)別反射波;后面兩個(gè)不易辨認(rèn)反射波。對(duì)每個(gè)初始模型都進(jìn)行了6次反演測(cè)試,包括1次常規(guī)反演和5次OT反演?？偣灿肧EAM Ⅱ Foothill模型進(jìn)行了24次反演實(shí)驗(yàn)。

圖12 實(shí)際速度模型

圖13 觀測(cè)數(shù)據(jù)

圖14 初始模型

圖15 根據(jù)圖14計(jì)算出來(lái)的合成地震數(shù)據(jù)

第1個(gè)初始模型最接近真實(shí)模型,且其初始合成數(shù)據(jù)已經(jīng)包含了大量反射事件。圖16顯示了使用第1個(gè)初始模型的反演結(jié)果。所有6個(gè)反演都給出了正確的結(jié)構(gòu),且反演結(jié)果基本相當(dāng)。常規(guī)全波形反演中的L2目標(biāo)函數(shù)的吸引域(attraction basin)較小,從另一個(gè)角度來(lái)說(shuō),L2的分辨率較高。圖16的反演結(jié)果顯示了OT函數(shù)也同樣具有較高的分辨率。與此相比,若是采用波動(dòng)方程層析(wave equation tomography),通常需要用L2接著往下反演;而OT則具有一步到位的潛力。

第2個(gè)初始模型比第1個(gè)更為光滑,大部分清晰的反射結(jié)構(gòu)已經(jīng)被平滑掉。其合成數(shù)據(jù)僅能識(shí)別一個(gè)大的反射事件。初至波則與第1個(gè)模型的初至波較為接近。原因?yàn)?①初至波穿透較淺;②初至波沿著淺地表傳播,其到時(shí)本身就是模型結(jié)構(gòu)的加權(quán)平均,所以對(duì)模型光滑不會(huì)導(dǎo)致初至太大的變化。圖17 展示的是6個(gè)目標(biāo)函數(shù)對(duì)應(yīng)的反演結(jié)果。常規(guī)反演沒(méi)有完全重構(gòu)速度模型,其最深處能識(shí)別的結(jié)構(gòu)在中間2km深處。與真實(shí)模型相比,其淺部的速度缺乏背景場(chǎng)的更新,更接近于偏移模式。在早期的全波形反演發(fā)展中,地球物理學(xué)家已經(jīng)意識(shí)到反射波對(duì)于中間波長(zhǎng)的速度模型是缺乏反應(yīng)的(intermediate-wavelength gap),所以用L2函數(shù)會(huì)導(dǎo)致全波形反演向著偏移的模式演化,并陷入局部極值。圖17所示結(jié)果實(shí)際上揭示了OT函數(shù)對(duì)于走時(shí)信息的提取能力。OT函數(shù)是基于模式匹配的,模式上的匹配是隱含著運(yùn)動(dòng)學(xué)信息(kinematic),如圖17和18所示。一維OT提取的是走時(shí)信息,高維OT除了提取走時(shí)信息外,還提取了事件沿著空間軸方向的移動(dòng)信息。這些運(yùn)動(dòng)學(xué)信息的提取對(duì)于模型低波數(shù)成分的恢復(fù)起了決定性作用。從圖17的5個(gè)OT結(jié)果上看,3個(gè)二維OT函數(shù)的反演結(jié)果略微優(yōu)于兩個(gè)一維OT函數(shù)的結(jié)果。在模型的左下方,兩個(gè)一維OT的結(jié)果未能清晰地重構(gòu)出低速帶。需要說(shuō)明的是,一維GS-OT的用戶輸入?yún)?shù)不大容易選取最優(yōu)輸入?yún)?shù),在特定模型上累計(jì)的經(jīng)驗(yàn)推廣到另一個(gè)完全不同的模型上往往并不會(huì)得到同樣的結(jié)果。從原理上講,GS-OT是具有良好的凸性的,比W2更優(yōu)越,圖17e中展示的反演結(jié)果或許還有進(jìn)一步提高的可能。

在作者積累的全波形反演經(jīng)驗(yàn)中,模型陷入局部極值后可能會(huì)有一些難以解釋的現(xiàn)象。圖18展示的是用第3個(gè)初始模型得到的6個(gè)反演結(jié)果。相比于圖17中的常規(guī)反演,圖18中的常規(guī)反演得到的模型更新更深。雖然兩者都是因全波形反演的偏移模式而陷入局部極值,但無(wú)法解釋為何圖18a會(huì)有更深的模型更新。與圖17類似,圖18中同樣是3個(gè)二維OT函數(shù)的反演結(jié)果略微優(yōu)于兩個(gè)一維OT函數(shù)的結(jié)果。并且GS-1D的反演結(jié)果確實(shí)比W2更好。

圖18 以圖14c作為初始模型得到的反演結(jié)果

圖19展示的是初始模型為速度隨深度線性增加模型對(duì)應(yīng)的6個(gè)反演結(jié)果,對(duì)于所有這些函數(shù),線性增加的模型可能離真實(shí)的模型太遠(yuǎn)。W2的反演可能是因?yàn)樵谙萑刖植繕O值的情形下步長(zhǎng)不合適而出現(xiàn)了特別大和特別小的速度結(jié)構(gòu)。其它5個(gè)反演結(jié)果應(yīng)當(dāng)是全波形反演的偏移模式,背景模型難以更新。值得注意的是,KR-OT對(duì)應(yīng)的反演結(jié)果模型迭代得最深,這可能跟KR-OT獨(dú)特的平衡能力相關(guān)。KR-OT的輸入信號(hào)即是數(shù)據(jù)殘差,在經(jīng)過(guò)KR-OT平衡后,來(lái)自于深部的反射信號(hào)得到了加強(qiáng),所以導(dǎo)致較深的模型迭代。

圖19 以圖14d作為初始模型得到的反演結(jié)果

3 討論

圖11中的常規(guī)全波形反演陷入局部極值,即使沒(méi)有真實(shí)模型參考也能識(shí)別是局部極值。但在實(shí)際的全波形反演中,應(yīng)用常規(guī)全波形反演后遇到的情形并不像圖11那么清晰。關(guān)于實(shí)際數(shù)據(jù)中出現(xiàn)的局部極值現(xiàn)象,一個(gè)更為準(zhǔn)確的推斷是,在部分區(qū)域模型恢復(fù)得很好,部分區(qū)域與真正的結(jié)構(gòu)差異較大。所以導(dǎo)致了全波形反演與某些測(cè)井資料對(duì)應(yīng)得很好,而與其它測(cè)井信息則差異較大。也正是因?yàn)檫@個(gè)原因,才凸顯出凸性更好的目標(biāo)函數(shù)的重要性,這樣即可在絕大部分區(qū)域都能成功避免局部極值。從目標(biāo)函數(shù)的研究歷史看,從最初的點(diǎn)對(duì)點(diǎn)的比較,逐步發(fā)展為道與道的比較?；诘辣容^的目標(biāo)函數(shù)包括:包絡(luò)函數(shù)、互相關(guān)、反卷積和互相關(guān)走時(shí)差等一系列函數(shù)。最優(yōu)化傳輸是基于道的整個(gè)數(shù)據(jù)體的模式比較,支持整個(gè)數(shù)據(jù)體之間的比較。從研究歷史看,基于數(shù)據(jù)整體的比較是一個(gè)發(fā)展方向,其優(yōu)勢(shì)在于可以在全局比較中正確匹配地震事件的偏差,是全局優(yōu)化策略。

最優(yōu)化傳輸函數(shù)的凸性取決于:①范數(shù)。通常會(huì)選擇1或2范數(shù)(p-Wasserstein函數(shù))。通?？梢哉J(rèn)為2范數(shù)比1范數(shù)更凸。此外,在大多數(shù)情況下,2范數(shù)更容易操作。相比之下,1范數(shù)在零點(diǎn)附近不是平滑的,比如,|x|在零附近的導(dǎo)數(shù)就是間斷的。不恰當(dāng)?shù)貙?范數(shù)合并到全波形反演中可能會(huì)導(dǎo)致伴隨源在某些點(diǎn)或者區(qū)域的不連續(xù)性。KR-OT函數(shù)雖然是基于1范數(shù),但并沒(méi)有出現(xiàn)這種不連續(xù)的現(xiàn)象,因?yàn)镵R-OT是在對(duì)偶空間中表示的,1范數(shù)被轉(zhuǎn)移到了約束方程中。②預(yù)處理。對(duì)于數(shù)據(jù)預(yù)處理需要謹(jǐn)慎對(duì)待。比如W2常用到的移動(dòng)加規(guī)則化,其移動(dòng)量不宜太大?；诜蔷€性變換的預(yù)處理有很多,在作者的經(jīng)驗(yàn)中,tanh表現(xiàn)不錯(cuò),而且用戶輸入?yún)?shù)容易調(diào)節(jié)。③時(shí)間空間坐標(biāo)是否引入最優(yōu)化傳輸目標(biāo)函數(shù),圖空間最優(yōu)化傳輸是將時(shí)間坐標(biāo)和空間坐標(biāo)顯式地引入,其結(jié)果也變得更凸,同時(shí),避免了負(fù)信號(hào)的處理。

具有全局優(yōu)化潛力的最優(yōu)化傳輸函數(shù)的計(jì)算量比常規(guī)函數(shù)要大。在實(shí)用化方面,一維的最優(yōu)化傳輸函數(shù)的計(jì)算量很小,可以應(yīng)用在二維和三維全波形反演中。若將三維數(shù)據(jù)分割成二維的矩陣,那么,本文提到的3個(gè)二維最優(yōu)化傳輸函數(shù)同樣可以應(yīng)用,計(jì)算代價(jià)在可接受范圍內(nèi)。若是直接處理三維數(shù)據(jù)體,計(jì)算量較大,可能需要一些加速策略。從目標(biāo)函數(shù)的性能來(lái)看,沒(méi)有哪一個(gè)函數(shù)可以處理所有的情形,但是有些函數(shù)會(huì)在大部分地質(zhì)環(huán)境中表現(xiàn)良好。使用SEAM Ⅱ Foothill模型的全波形反演實(shí)驗(yàn)證明,總體趨勢(shì)是2D最優(yōu)化傳輸函數(shù)優(yōu)于1D版本。BFM功能優(yōu)于W2功能;GS-2D技術(shù)優(yōu)于GS-1D技術(shù)。從理論上說(shuō),二維函數(shù)能提取空間軸方向的信息,比一維更具潛力。這與作者在其它各種模型上測(cè)試得到的認(rèn)識(shí)是一致的。綜合作者在各種模型上的測(cè)試經(jīng)驗(yàn)看,W2在數(shù)據(jù)模式較為簡(jiǎn)單的情形下表現(xiàn)良好,可能是因?yàn)轭A(yù)處理影響到了函數(shù)性能。KR-OT因?yàn)檩斎霐?shù)據(jù)是數(shù)據(jù)誤差,導(dǎo)致地震事件從一開(kāi)始就混雜在一起,所以合成數(shù)據(jù)跟真實(shí)數(shù)據(jù)不能相差太遠(yuǎn)。但是其具有獨(dú)特且優(yōu)越的振幅平衡能力。在陸地?cái)?shù)據(jù)的彈性全波形反演中KR-OT是首選。從原理上說(shuō),KR-OT比L2更凸,但會(huì)比GS-1D略遜一籌。GS-OT(包括GS-1D和GS-2D)的優(yōu)勢(shì)在于將坐標(biāo)引入計(jì)算代價(jià)中,從根本上避免了數(shù)據(jù)的預(yù)處理。但是,凸性良好的函數(shù)具有共同的小瑕疵,即在反演模型接近真實(shí)模型時(shí),反映在數(shù)據(jù)上的偏差可能不會(huì)太大,所以反演收斂速度會(huì)很慢?；诖?在GS-OT反演后,可以酌情在后續(xù)加入常規(guī)反演。BFM的性能原則上比KR-OT更凸,但遜于GS-2D。BFM是基于2范數(shù)的最優(yōu)化傳輸,而KR-OT是基于1范數(shù)的,2范數(shù)比1范數(shù)更凸。GS-2D也是基于2范數(shù),同時(shí),避免了數(shù)據(jù)預(yù)處理,應(yīng)當(dāng)比BFM更凸。

4 結(jié)論

本文綜合對(duì)比了目前學(xué)術(shù)界和工業(yè)界新提出的最優(yōu)化傳輸目標(biāo)函數(shù)。常規(guī)的最小二乘目標(biāo)函數(shù)是點(diǎn)對(duì)點(diǎn)的比較方式,對(duì)于勘探地震數(shù)據(jù)地震道中的時(shí)間相關(guān)性和地震道之間的空間相關(guān)性沒(méi)有顯式地提取。當(dāng)初始模型足夠準(zhǔn)確,合成地震事件與真實(shí)地震事件偏差在半個(gè)波長(zhǎng)以內(nèi)的時(shí)候,常規(guī)的最小二乘目標(biāo)函數(shù)是可以反演出足夠接近真實(shí)的地下速度模型。然而在實(shí)際數(shù)據(jù)處理中,相對(duì)于采集中設(shè)定的低頻范圍,初始模型常常是不夠準(zhǔn)確的。這也是學(xué)術(shù)界和工業(yè)界尋求其它目標(biāo)函數(shù)的動(dòng)力。從目標(biāo)函數(shù)的研究歷史看,的確是朝著高維方向發(fā)展。在目前的研究中,最優(yōu)化傳輸通常采用1范數(shù)和2范數(shù);對(duì)于給定的范數(shù),又可以分為單道處理和高維處理。綜合我們的測(cè)試結(jié)果看,高維的最優(yōu)化傳輸函數(shù)優(yōu)于一維的版本。在實(shí)際應(yīng)用中,可以綜合考慮計(jì)算成本和效果。對(duì)于工區(qū)不太復(fù)雜的情形,可以采用W2,因?yàn)槠浯嬖诮馕鼋?計(jì)算速度最快。在復(fù)雜地區(qū),可以采用BFM,KR-OT,GS-1D或者GS-2D;在這4種最優(yōu)化傳輸函數(shù)中,以GS-1D計(jì)算成本最小。其中,只有KR-OT是基于1范數(shù),其凸性可能遜于BFM;但是KR-OT獨(dú)有的振幅平衡能力使得KR-OT在處理存在強(qiáng)振幅差異地震數(shù)據(jù)的時(shí)候具有特殊的優(yōu)勢(shì)。