本發(fā)明涉及一種地震勘探的地下成像方法，尤其是地震干涉法被動源低頻重構(gòu)地震超越干涉法主動源數(shù)據(jù)、分炮并行直達(dá)波估計和數(shù)據(jù)驅(qū)動面向復(fù)雜高速屏蔽地下構(gòu)造的高精度marchenko成像。
背景技術(shù)：
：：在地震成像反演過程中，低頻信息對于地下介質(zhì)的構(gòu)造成像尤為重要，它直接影響地下介質(zhì)精細(xì)結(jié)構(gòu)的反演效果。低頻信號具有較強(qiáng)的穿透高速屏蔽層的能力，可以利用低頻信號來降低高速屏蔽層對地震信號的散射和屏蔽作用，提高深層信號的信噪比。實際應(yīng)用中，由于被動源數(shù)據(jù)具有固有的低頻優(yōu)勢，常常通過地震干涉法對被動源數(shù)據(jù)進(jìn)行處理，得到有用的低頻信息。地震干涉理論由claerbout于1968年提出，calvert于2004基于地震干涉理論提出“虛震源法”，可以將震源重置到檢波點位置，該方法在勘探地震數(shù)據(jù)上得到了廣泛應(yīng)用。wapenaar于2012年用積分理論和互易定理證明了地震干涉技術(shù)對不同的非衰減介質(zhì)和不同的震源條件下均成立，并用理論和實際數(shù)據(jù)進(jìn)行了驗證。作為地震干涉法的拓展，地震超越干涉法最早由keeswapenaar和filippobroggini于2012年提出，此方法實現(xiàn)了在無需知道地下介質(zhì)信息的情況下在非檢波器位置構(gòu)建虛擬震源的新方法，即虛擬震源的位置可以分布在地下介質(zhì)的任意位置。地震超越干涉法利用地表接收到的主動源反射響應(yīng)和在不精確的平滑模型上估計的直達(dá)波來重構(gòu)包含一次波和多次波的格林函數(shù)場。利用重構(gòu)得到的上傳和下傳格林函數(shù)進(jìn)行marchenko成像，能夠在無需了解上覆介質(zhì)構(gòu)造的情況準(zhǔn)確對目標(biāo)區(qū)域的介質(zhì)進(jìn)行精確成像?；诘皖l缺失的地震數(shù)據(jù)重構(gòu)地下介質(zhì)處的格林函數(shù)會導(dǎo)致在最終的marchenko成像結(jié)果中出現(xiàn)由于低頻缺失旁瓣導(dǎo)致的地下反射界面不能聚焦的問題，這會對后續(xù)地震勘探的解釋造成影響。實際中對地震數(shù)據(jù)低頻缺失的重構(gòu)辦法有兩種，一種是從數(shù)學(xué)計算的方式得到低頻信息：han于2012年和2014年基于壓縮感知和稀疏反演方法重構(gòu)數(shù)據(jù)低頻信息；huang于2015年用彈性波包絡(luò)法恢復(fù)了地下介質(zhì)縱橫波速度的長波長成分。此類低頻重構(gòu)方法在實際應(yīng)用中可能在地震記錄上產(chǎn)生虛假信息，最終影響成像結(jié)果。另外一種是對被動源數(shù)據(jù)進(jìn)行地震干涉法處理，然后對主動源進(jìn)行低頻重構(gòu)：zhang于2015年從能量角度(振幅)匹配主動源和被動源數(shù)據(jù)，但即使經(jīng)過去噪之后的高信噪比低頻數(shù)據(jù)中依然存在噪聲，這種局限性也會對最終的成像結(jié)果產(chǎn)生影響。另外，此種被動源低頻補(bǔ)充主動源數(shù)據(jù)方法是基于人工增益系數(shù)來完成的，主觀性較強(qiáng)，且不能基于數(shù)據(jù)驅(qū)動。marchenko成像算法依賴于目標(biāo)成像區(qū)域的每一個點的格林函數(shù)，很小的網(wǎng)格大小會造成巨大的計算量。傳統(tǒng)的單炮串行計算遠(yuǎn)遠(yuǎn)不能滿足計算效率的要求。技術(shù)實現(xiàn)要素：：本發(fā)明的目的針對上述現(xiàn)有技術(shù)的不足，提供一種提高低頻重構(gòu)效果和成像效率的并行面向目標(biāo)數(shù)據(jù)驅(qū)動marchenko成像方法。本發(fā)明的目的是通過以下技術(shù)方案實現(xiàn)的：低頻重構(gòu)并行面向目標(biāo)數(shù)據(jù)驅(qū)動marchenko成像方法的核心是基于頻率優(yōu)勢的地震干涉法被動源數(shù)據(jù)對地震超越干涉法主動源數(shù)據(jù)低頻重構(gòu)、基于pbs批量作業(yè)管理系統(tǒng)的并行地震超越干涉法格林函數(shù)重構(gòu)及面向復(fù)雜地下構(gòu)造的marchenko成像。對復(fù)雜地下構(gòu)造的精確marchenko成像結(jié)果需要寬頻帶的地震數(shù)據(jù)作為輸入，首先基于地震干涉法被動源數(shù)據(jù)的低頻優(yōu)勢對地震超越干涉法進(jìn)行低頻重構(gòu)，然后基于歸一化加和之后根據(jù)主動源數(shù)據(jù)最大振幅進(jìn)行能量恢復(fù)，得到寬頻帶的地震數(shù)據(jù)。迭代求取聚焦函數(shù)、完全格林函數(shù)、上傳格林函數(shù)和下傳格林函數(shù)的計算量巨大，使用pbs批量作業(yè)管理系統(tǒng)將計算平臺分為一個主處理器和若干從處理器。主處理器監(jiān)控整個平滑模型上估計直達(dá)波分炮正演模擬階段，由主處理器將正演計算分配到各個從處理器上進(jìn)行，計算完成之后再由主處理器收集完全格林函數(shù)重構(gòu)結(jié)果、完全格林函數(shù)分解結(jié)果，然后提交基于互相關(guān)marchenko成像的處理作業(yè)到從處理器進(jìn)行計算，此時主處理器監(jiān)控整個基于上傳和下傳格林函數(shù)計算虛擬震源位置反射系數(shù)階段，計算完成之后，主處理器收集得到的目標(biāo)區(qū)域每一個地下虛擬震源位置的反射系數(shù)，再根據(jù)地下虛擬震源的位置對反射系數(shù)按所處的觀測系統(tǒng)位置依次排列，得到最終的marchenko成像結(jié)果。這樣既解決了地下復(fù)雜地質(zhì)構(gòu)造的高精度成像問題，又基于分炮并行正演充分發(fā)揮了當(dāng)前計算平臺的優(yōu)秀性能，顯著縮短了marchenko成像的計算時間。低頻重構(gòu)并行marchenko成像方法，包括以下步驟：a、在ubuntu14.04系統(tǒng)平臺上安裝torque2.4.16，構(gòu)建pbs批量作業(yè)管理系統(tǒng)；b、對實際采集的被動源地震信號進(jìn)行地震干涉法處理，基于互相關(guān)法，重構(gòu)與主動源相同觀測系統(tǒng)的被動源地震數(shù)據(jù)；c、對被動源地震數(shù)據(jù)進(jìn)行多域迭代去噪，得到較高信噪比的被動源數(shù)據(jù)。補(bǔ)償去噪之后的被動源數(shù)據(jù)能量，找到被動源低頻端與主動源低頻端在頻譜上的交點，并確定焦點位置的頻率為優(yōu)勢頻率；d、結(jié)合低通濾波器分別對主動源、被動源在優(yōu)勢頻率處做低通濾波；e、設(shè)計匹配濾波器，計算匹配濾波算子，把低通濾波后的主動源低頻信息與低通濾波后的被動源低頻信息匹配濾波，得到匹配濾波后的主動源低頻信息；f、對匹配濾波后的主動源低頻信息和主動源分別做歸一化，然后對二者做加和，然后基于主動源歸一化之前的主動源數(shù)據(jù)最大振幅進(jìn)行能量恢復(fù)，得到合成后的寬頻帶地震數(shù)據(jù)；g、對地震數(shù)據(jù)預(yù)處理，滿足marchenko迭代求解格林函數(shù)的要求，這個步驟分兩種情況：a、數(shù)值模擬算例中，通過數(shù)值計算的方法來計算出寬頻帶地震數(shù)據(jù)中的直達(dá)波、自由表面多次波。b、在實際地震勘探中，自由表面多次波的消去可以通過srme方法達(dá)到，直達(dá)波的計算可以通過淺層調(diào)查進(jìn)行計算得到；h、在不精確的模型上進(jìn)行正演得到地下介質(zhì)某點到地表布設(shè)的檢波器之間的估計的直達(dá)波響應(yīng)。經(jīng)過低頻補(bǔ)償后的寬頻帶地震數(shù)據(jù)可以對地下復(fù)雜地質(zhì)區(qū)域進(jìn)行成像。選定地下介質(zhì)由于類似于高阻抗地層的屏蔽作用導(dǎo)致反射地震能量較弱、常規(guī)的地震成像方法不夠精確的地區(qū)作為成像區(qū)域。假設(shè)成像區(qū)域有m個虛擬震源位置，那么就需要進(jìn)行m次分炮正演。假設(shè)計算節(jié)點上有n個從處理器參與計算，那么把并行計算的作業(yè)提交到pbs作業(yè)系統(tǒng)后n個從處理器可以進(jìn)行并行計算，當(dāng)一炮的正演結(jié)束后，從處理器結(jié)束計算，此時，pbs主處理器為這個從處理器按順序分配下一個正演計算作業(yè)，直到m次正演全部結(jié)果；i、marchenko成像需要求解成像區(qū)域的所有成像點的上傳格林函數(shù)和下傳格林函數(shù)。成像區(qū)域m個虛擬震源的位置對應(yīng)著m次迭代求解格林函數(shù)。pbs主處理器將經(jīng)過預(yù)處理之后的寬頻帶地震數(shù)據(jù)和h步里正演得到估計的直達(dá)波分配到n個從處理當(dāng)中。對于每一個從處理器的單次marchenko迭代求解來說，首先迭代求解三維marchenko方程得到聚焦函數(shù)，再由完全格林函數(shù)與聚焦函數(shù)2、寬頻帶數(shù)據(jù)之間的關(guān)系，求得完全格林函數(shù)，之后依據(jù)上傳、下傳格林函數(shù)與聚焦函數(shù)1之間的關(guān)系，求得分解后的上傳格林函數(shù)和下傳格林函數(shù)。當(dāng)一次marchenko迭代求解完全格林函數(shù)、上傳格林函數(shù)和下傳格林函數(shù)結(jié)束，主處理器將繼續(xù)分配后續(xù)作業(yè)給從處理器處理，直到全部的作業(yè)結(jié)束；j、將分解之后得到的上傳格林函數(shù)和下傳格林函數(shù)用來進(jìn)行marchenko成像。成像區(qū)域m個虛擬震源位置對應(yīng)著m次單個虛擬震源位置的上傳格林函數(shù)與下傳格林函數(shù)的互相關(guān)函數(shù)計算以及該位置的反射系數(shù)求取。pbs主處理器首先將第1,2,3，…，n個計算作業(yè)依次分配到第1,2,3，…，n個從處理器進(jìn)行計算。當(dāng)單個虛擬震源位置的反射系數(shù)計算得到之后，pbs主處理器收集記錄該虛擬震源位置的反射系數(shù)并記錄該虛擬震源在觀測系統(tǒng)中的坐標(biāo)，然后為空閑的從處理器繼續(xù)分配作業(yè)，直到全部的作業(yè)計算完成；k、當(dāng)系統(tǒng)求得成像區(qū)域內(nèi)全部虛擬震源位置的反射系數(shù)后。根據(jù)觀測系統(tǒng)，依據(jù)成像區(qū)域中的網(wǎng)格點坐標(biāo)的先后順序，串行依次把反射系數(shù)排列到網(wǎng)格點的位置，最終得到基于低頻重構(gòu)后的面向目標(biāo)區(qū)域的數(shù)據(jù)驅(qū)動marchenko成像結(jié)果。有益效果：本發(fā)明成功地將基于頻率優(yōu)勢的被動源低頻信息重構(gòu)主動源低頻信息技術(shù)應(yīng)用到了marchenko成像中，提出的基于被動源低頻優(yōu)勢的主動源低頻重構(gòu)方法基于匹配濾波器能夠最大限度地防止被動源低頻噪聲對重構(gòu)出的寬頻帶地震數(shù)據(jù)的影響，同時基于歸一化能量恢復(fù)方法既恢復(fù)得到了寬頻帶地震數(shù)據(jù)，又完整地保留了主動源數(shù)據(jù)的高頻部分，這比人為地控制與調(diào)整低頻信號在寬頻帶地震數(shù)據(jù)中的權(quán)重更符合現(xiàn)實，且本方法基于數(shù)據(jù)驅(qū)動。提出的基于重構(gòu)后的寬頻帶地震數(shù)據(jù)分炮并行迭代求解聚焦函數(shù)、完全格林函數(shù)、上傳格林函數(shù)、下傳格林函數(shù)，節(jié)約了大量的時間，顯著提高了marchenko成像的計算效率，最終在寬頻帶地震數(shù)據(jù)marchenko成像中地下反射界面更加聚焦，相對傳統(tǒng)的低頻補(bǔ)償方法，本方面沒有引入虛假的反射界面。得到了更好的成像效果。本發(fā)明主要解決了以下問題：1、地震干涉法和地震超越干涉法基于不同的假設(shè)，本發(fā)明分別從數(shù)值模擬的角度和實際生產(chǎn)的角度，給出了把此兩種數(shù)據(jù)結(jié)合的方法。為被動源數(shù)據(jù)在marchenko成像中的應(yīng)用奠定了基礎(chǔ)。2、常規(guī)采集的主動源地震數(shù)據(jù)中含有的低頻信號往往弱于其高頻信號，寬頻帶的地震數(shù)據(jù)需要基于低頻重構(gòu)得到。相對于純數(shù)值低頻延拓，本發(fā)明所用的低頻信號來源于真實的地震數(shù)據(jù)，所以補(bǔ)償?shù)玫降膶掝l帶地震數(shù)據(jù)更加可靠；相對于傳統(tǒng)的基于能量的補(bǔ)償方法會在最終的成像結(jié)果中產(chǎn)生虛假反射界面，本發(fā)明提出的頻率優(yōu)勢地重構(gòu)方法能夠在完整保留主動源高頻信號的基礎(chǔ)上，很好地重構(gòu)了主動源的低頻能量，并且被動源的噪聲信號不會對最終的成像結(jié)果造成影響。3、傳統(tǒng)的地震方法對含有高阻抗上覆地層的區(qū)域往往由于強(qiáng)散射、高速屏蔽作用導(dǎo)致照明不足、成像效果較差。本發(fā)明同時利用了低頻信息較強(qiáng)的高速屏蔽層穿透能力和marchenko成像算法對層間多次波的消去作用，獲得了高精度的地下構(gòu)造成像結(jié)果。4、marchenko成像算法依賴于目標(biāo)成像區(qū)域的每一個點的格林函數(shù)，很小的網(wǎng)格大小會造成巨大的計算量，通過pbs批量作業(yè)管理系統(tǒng)把計算作業(yè)分配到從處理器上進(jìn)行計算，可以得到線性的加速比。附圖說明圖1低頻重構(gòu)并行marchenko成像方法流程圖。圖2模型圖(左)速度模型，(右)密度模型圖3不精確的模型圖(左)平滑后的速度模型，(右)平滑后的密度模型圖4目標(biāo)成像區(qū)域示意圖(左)目標(biāo)成像區(qū)域，(右)放大后的目標(biāo)成像區(qū)域圖5被動源互相關(guān)地震干涉法處理之后的單炮記錄圖6主動源單炮記錄和低頻重構(gòu)后的寬頻帶單炮記錄(左)主動源單炮記錄，(右)低頻重構(gòu)后的寬頻帶單炮記錄。圖7低頻重構(gòu)后抽道波形和全頻帶正演抽道波形圖8低頻缺失單炮頻譜、低頻重構(gòu)后寬頻帶單炮頻譜和全頻帶正演單炮頻譜圖9主動源數(shù)據(jù)重構(gòu)的格林函數(shù)圖10主動源低頻重構(gòu)后重構(gòu)的格林函數(shù)圖11參考格林函數(shù)圖12主動源數(shù)據(jù)重構(gòu)的格林函數(shù)、主動源低頻重構(gòu)后重構(gòu)的格林函數(shù)與參考格林函數(shù)抽道對比圖13marchenko成像結(jié)果(左)低頻缺失主動源地震數(shù)據(jù)的marchenko成像結(jié)果，(右)低頻重構(gòu)之后的寬頻帶地震數(shù)據(jù)的marchenko成像結(jié)果。圖14相同區(qū)域逆時偏移成像結(jié)果具體實施方式下面結(jié)合附圖和實例對本發(fā)明進(jìn)一步的詳細(xì)說明低頻重構(gòu)并行marchenko成像方法，包括以下步驟：a、程序是基于c語言編寫完成，為了進(jìn)行分炮并行，在ubuntu14.04系統(tǒng)平臺上安裝torque2.4.16，構(gòu)建pbs批量作業(yè)管理系統(tǒng)；b、對實際采集的被動源地震信號進(jìn)行地震干涉法處理，基于下列互相關(guān)法，重構(gòu)與主動源相同觀測系統(tǒng)的被動源地震數(shù)據(jù)，為了增強(qiáng)重構(gòu)之后地震數(shù)據(jù)的信噪比，需要長時間地接收被動源信號。其中分別為a、b位置處的觀測波場，表示xb為震源、xa為檢波器的格林函數(shù)，為功率譜，代表實部，c和ρ分別代表傳播速度和密度；c、匹配濾波需要輸入高信噪比的地震數(shù)據(jù)，需要對被動源地震數(shù)據(jù)進(jìn)行多域迭代去噪得到較高信噪比的被動源數(shù)據(jù)。補(bǔ)償去噪之后的被動源數(shù)據(jù)能量，找到被動源低頻端與主動源低頻端在頻譜上的交點，并確定焦點位置的頻率為優(yōu)勢頻率fcut；d、結(jié)合低通濾波器對主動源做低通濾波得到低頻濾波后的主動源地震數(shù)據(jù)、對被動源在優(yōu)勢頻率處做低通濾波得到低頻濾波后的被動源地震數(shù)據(jù)：e、設(shè)計匹配濾波器，計算匹配濾波算子m(t)，把低通濾波后的主動源低頻信息與低通濾波后的被動源低頻信息匹配濾波，得到匹配濾波后的主動源低頻信息：f、對匹配濾波后的主動源低頻信息和主動源分別做歸一化，然后對二者做加和，然后基于主動源歸一化之前的主動源數(shù)據(jù)最大振幅進(jìn)行能量恢復(fù)，得到合成后的寬頻帶地震數(shù)據(jù):g、對地震數(shù)據(jù)預(yù)處理，滿足marchenko迭代求解格林函數(shù)的要求，這個步驟分兩種情況：a、在數(shù)值模擬算例中，通過數(shù)值計算的方法來計算出寬頻帶地震數(shù)據(jù)中的直達(dá)波、自由表面多次波:假設(shè)f1(t)是地表存在吸收邊界條件(如pml)的反射響應(yīng)，其包含了一次波、直達(dá)波和層間多次波。f2(t)是在均勻介質(zhì)數(shù)值模擬得到的直達(dá)波響應(yīng)。那么自由表面多次波f3(t)可以由下式得出：f3(t)＝f(t)-f1(t)此時符合marchenko成像的反射響應(yīng)r0可以表示為：r0＝f(t)-f3(t)-f2(t)b、在實際地震勘探中，自由表面多次波的消去可以通過srme方法達(dá)到，直達(dá)波的計算可以通過淺層調(diào)查進(jìn)行計算得到；h、在不精確的模型上(附圖3)進(jìn)行正演得到地下介質(zhì)某點到地表布設(shè)的檢波器之間的估計的直達(dá)波響應(yīng)。值得注意的是，此時的模型不需要精確的模型進(jìn)行直達(dá)波估計。經(jīng)過低頻補(bǔ)償后的寬頻帶地震數(shù)據(jù)可以對地下復(fù)雜地質(zhì)區(qū)域進(jìn)行成像。選定地下介質(zhì)由于類似于高阻抗地層的屏蔽作用導(dǎo)致反射地震能量較弱、常規(guī)的地震成像方法不夠精確的地區(qū)作為成像區(qū)域(附圖4)。假設(shè)成像區(qū)域有m個虛擬震源位置，那么就需要進(jìn)行m次分炮正演。假設(shè)計算節(jié)點上有n個從處理器參與計算，那么把并行計算的作業(yè)提交到pbs作業(yè)系統(tǒng)后n個從處理器可以進(jìn)行并行計算，當(dāng)一炮的正演結(jié)束后，從處理器結(jié)束計算，此時，pbs主處理器為這個從處理器按順序分配下一個正演計算作業(yè)，直到m次正演全部結(jié)果；i、marchenko成像需要求解成像區(qū)域的所有成像點的上傳格林函數(shù)和下傳格林函數(shù)。上傳格林函數(shù)g-、下傳格林函數(shù)g+與聚焦函數(shù)及反射地震數(shù)據(jù)的關(guān)系為其中r為反射地震數(shù)據(jù)，f1+和f1-分別為上傳和下傳聚焦函數(shù)，t代表時間。依據(jù)上述表達(dá)式，已知聚焦函數(shù)，可以求取上傳格林函數(shù)和下傳格林函數(shù)。聚焦函數(shù)需要通過迭代求取三維marchenko方程得到。一般迭代次數(shù)設(shè)為4-15次。成像區(qū)域m個虛擬震源的位置對應(yīng)著m次迭代求解格林函數(shù)及其分解。pbs主處理器將經(jīng)過預(yù)處理之后的寬頻帶地震數(shù)據(jù)和h步里正演得到估計的直達(dá)波分配到n個從處理當(dāng)中。對于每一個從處理器的單次marchenko迭代求解來說，首先迭代求解三維marchenko方程得到聚焦函數(shù)，再由上傳、下傳格林函數(shù)與聚焦函數(shù)f1+、f1-之間的關(guān)系，求得上傳格林函數(shù)和下傳格林函數(shù)。當(dāng)一次marchenko迭代求解上傳格林函數(shù)和下傳格林函數(shù)結(jié)束，主處理器將繼續(xù)分配后續(xù)作業(yè)給從處理器處理，直到全部的作業(yè)結(jié)束；j、將分解之后得到的上傳格林函數(shù)g-和下傳格林函數(shù)g+用來進(jìn)行marchenko成像。成像區(qū)域m個虛擬震源位置對應(yīng)著m次單個虛擬震源位置的上傳格林函數(shù)g-與下傳格林函數(shù)g+的互相關(guān)函數(shù)計算以及該位置的反射系數(shù)求取?；ハ嚓P(guān)函數(shù)表示為：其中，x″0為面上檢波器。在頻率域，互相關(guān)成像條件是互相關(guān)函數(shù)在所有頻率的積分。pbs主處理器首先將第1,2,3，…，n個計算作業(yè)依次分配到第1,2,3，…，n個從處理器進(jìn)行計算。當(dāng)單個虛擬震源位置的反射系數(shù)計算得到之后，pbs主處理器收集記錄該虛擬震源位置的反射系數(shù)并記錄該虛擬震源在觀測系統(tǒng)中的坐標(biāo)，然后為空閑的從處理器繼續(xù)分配作業(yè)，直到全部的作業(yè)計算完成；k、當(dāng)系統(tǒng)求得成像區(qū)域內(nèi)全部虛擬震源位置的反射系數(shù)后。根據(jù)觀測系統(tǒng)，依據(jù)成像區(qū)域中的網(wǎng)格點坐標(biāo)的先后順序，串行依次把反射系數(shù)排列到網(wǎng)格點的位置，最終得到基于低頻重構(gòu)后的面向目標(biāo)區(qū)域的數(shù)據(jù)驅(qū)動marchenko成像結(jié)果。實施例1根據(jù)地震勘探成像的計算要求，torque2.4.16在ubuntu14.04系統(tǒng)下進(jìn)行安裝，進(jìn)行pbs并行平臺的搭建。利用含有高阻抗高速屏蔽褶皺地層的速度和密度模型(附圖2)進(jìn)行測試。模型參數(shù)如下：模型網(wǎng)格大小為801×2401，網(wǎng)格距dz＝dx＝2.5m，橫向距離為6.0km，縱向深度為2.0km，模型中地震波速度范圍從1.9km/s到2.8km/s，地震檢波器安置在模型表面，共601個檢波器，且檢波器之間的間隔為10m，主動源正演模擬表面震源位置同檢波器位置，共有601個震源。檢波器采樣間隔為0.004s。為了滿足互相關(guān)法重構(gòu)的假設(shè)，主動源和被動源的地表均是一個自由表面。首先在模型上模擬被動源，為了得到高信噪比的被動源重構(gòu)結(jié)果，模擬地下4500個噪聲源隨機(jī)分布在橫向(-2900m,2900m)，縱向(1350m,1700m)范圍內(nèi)，檢波器接收1200s的被動源數(shù)據(jù)。通過常規(guī)互相關(guān)重構(gòu)進(jìn)行被動源重構(gòu)，地表(0m,0m)位置的虛擬震源的地震記錄如附圖5。然后進(jìn)行主動源模擬，模擬得到20hz以下低頻信號較弱的的主動源地震數(shù)據(jù)(附圖6左)。根據(jù)本發(fā)明提出的頻率優(yōu)勢低頻重構(gòu)方法，首先經(jīng)過對被動源的能量補(bǔ)償，確定20hz為優(yōu)勢頻率。結(jié)合低通濾波器，分別對主動源和被動源信號進(jìn)行20hz低通濾波，分別得到低頻濾波之后的主動源數(shù)據(jù)和被動源數(shù)據(jù)。將此時的主動源數(shù)據(jù)去匹配被動源數(shù)據(jù)：首先求解匹配濾波算子，然后將匹配濾波算子作用于低通濾波之后的主動源數(shù)據(jù)，得到重構(gòu)出的主動源低頻信息。然后結(jié)合歸一化方法，對地震數(shù)據(jù)進(jìn)行能量恢復(fù)，得到含有直達(dá)波、自由表面多次波、層間多次波的寬頻帶地震數(shù)據(jù)(附圖6右)，由附圖7可知，重構(gòu)后的寬頻帶地震數(shù)據(jù)與全頻帶正演得到的地震數(shù)據(jù)符合地很好。對主動源數(shù)據(jù)、重構(gòu)后的寬頻帶地震數(shù)據(jù)以及全頻帶正演數(shù)據(jù)進(jìn)行頻譜分析(附圖8)可知，原始主動源中的高頻成分得到了很好的保留，同時由低頻端的頻譜趨勢可知低頻信號得到了有效的重構(gòu)。地震超越干涉法重構(gòu)格林函數(shù)得輸入是不包含直達(dá)波、自由表面多次波的數(shù)據(jù)，需要對數(shù)據(jù)進(jìn)行預(yù)處理?；诒景l(fā)明提出的數(shù)值模擬預(yù)處理方法，對寬頻帶地震數(shù)據(jù)進(jìn)行預(yù)處理后作為地震超越干涉法的輸入。對直達(dá)波的估計在不精確的模型上進(jìn)行(附圖3)。然后在目標(biāo)成像區(qū)域(附圖4)基于并行的方式迭代求取完全格林函數(shù)、上傳格林函數(shù)、下傳格林函數(shù)。此時需要給定迭代次數(shù)，一般的迭代次數(shù)設(shè)定為4-15之間。(0,1000)處的主動源數(shù)據(jù)重構(gòu)的格林函數(shù)(附圖9)、主動源低頻補(bǔ)償后重構(gòu)的格林函數(shù)(附圖10)、參考格林函數(shù)(附圖11)如附圖所示。對此三種格林函數(shù)抽道對比如圖所示，可以看到低頻重構(gòu)之后的格林函數(shù)與參考格林函數(shù)更加符合。通過對比低頻缺失marchenko成像結(jié)果與低頻重構(gòu)marchenko成像結(jié)果，我們可以看到，低頻重構(gòu)的反射界面更加聚焦(附圖13右)，模型里的具有高速屏蔽的褶皺的界面得到了很好的成像，深層界面成像清晰。而在低頻缺失marchenko成像(附圖13左)中，由于輸入的反射響應(yīng)里低頻缺失，導(dǎo)致旁瓣效應(yīng)影響了marchenko成像效果。最終的結(jié)果驗證了本發(fā)明的有效性。由于最終的成像結(jié)果中沒有由低頻匹配濾波引入的虛假反射界面，所以本發(fā)明中所提出的方法具有抗噪性。同時，在相同區(qū)域與傳統(tǒng)的成像方法的成像結(jié)果進(jìn)行對比可知本方法所提出的方法在成像上的優(yōu)勢。表1性能測試環(huán)境表2單核及多核并行計算效率對比結(jié)果從處理器核心數(shù)計算耗時cpu使用率最大使用內(nèi)存(mb)124分29秒26.1％91.5mb214分53秒52.1％184mb47分10秒98.6％370mb利用基于表1的測試環(huán)境，我們分別對marchenko成像過程進(jìn)行了三組測試。分別計算11個地下虛擬震源位置的marchenko成像結(jié)果。測試結(jié)果見表2，由表2可以看出，在包含直達(dá)波估計、求解格林函數(shù)以及基于互相關(guān)函數(shù)的成像的過程，隨著核心數(shù)量的增加，在表1所示的測試環(huán)境下，獲得了接近于線性的加速比。圖1是整個低頻重構(gòu)并行marchenko成像流程圖，從流程圖可以看出低頻重構(gòu)并行marchenko方法共分為以下幾個步驟：(1)基于互相關(guān)法重構(gòu)與主動源相同觀測系統(tǒng)的被動源地震數(shù)據(jù)；(2)補(bǔ)償去噪之后的被動源數(shù)據(jù)能量，確定優(yōu)勢頻率f_cut；(3)結(jié)合低通濾波進(jìn)行主動源低頻與被動源低頻匹配濾波；(4)歸一化能量匹配合成寬頻帶地震數(shù)據(jù)；(5)進(jìn)行數(shù)據(jù)預(yù)處理，利用數(shù)值模擬方法消去自由表面多次波及直達(dá)波后將其以及估計的地下虛擬震源到地表檢波器間的直達(dá)波作為輸入，并基于marchenko方程迭代并行求解格林函數(shù)；(6)最后基于互相關(guān)成像條件，循環(huán)求取目標(biāo)區(qū)域的所有成像點位置的反射系數(shù)完成marchenko成像。當(dāng)前第1頁12