地球物理反演與層析成像

地球物理反演與層析成像結(jié)業(yè)論文及地震走時(shí)層析成像程序姓名：學(xué)號(hào)：班級(jí)：地球物理學(xué)*班［在此處鍵入文檔的摘要。摘要通常是對(duì)文檔內(nèi)容的簡(jiǎn)短總結(jié)。在此處鍵入文檔的摘要。摘要通常是對(duì)文檔內(nèi)容的簡(jiǎn)短總結(jié)。］概述地學(xué)層析成像是用醫(yī)學(xué)X射線CT的理論詳細(xì)調(diào)查地下物性參數(shù)分布狀況的物探技術(shù)。分為地震層析成像、電磁波層析成像和電阻率層析成像。地震層析成像就是用地震數(shù)據(jù)來反演地下結(jié)構(gòu)的物質(zhì)屬性，并逐層剖析繪制其圖像的技術(shù)。其主要目的是確定地球內(nèi)部的精細(xì)結(jié)構(gòu)和局部不均勻性。相對(duì)來說，地震層析成像較其他兩種方法應(yīng)用更加廣泛，這是因?yàn)榈卣鸩ǖ乃俣扰c巖石性質(zhì)有比較穩(wěn)定的相關(guān)性，地震波衰減程度比電磁波小，且電磁波速度快，不易測(cè)量。地震層析成像按研究區(qū)域的尺度可分為全球?qū)游龀上瘛^(qū)域?qū)游龀上?、局部層析成像：按所用資料的來源可分為天然地震層析成像（大尺度深部橫向不均勻性研究）、人工地震測(cè)深（主要研究淺部界面分布）。按所依據(jù)的理論基礎(chǔ)一般分為基于射線方程的層析成像和基于波動(dòng)方程的層析成像。前者按射線追蹤時(shí)所用的地震波資料的不同又可分為體波（反射波、折射波）和面波層析成像：按反演的物性參數(shù)區(qū)分，可分為利用地震波走時(shí)反演地震波速度的波速層析成像以及利用地震波振幅衰減反演地震波衰減系數(shù)的層析成像。基于射線理論，地震波走時(shí)層析成像方法由于走時(shí)具有較高信噪比、無論是柱面波還是球面波走時(shí)的規(guī)律都相同等優(yōu)點(diǎn)，相對(duì)來說發(fā)展較早，技術(shù)方法比較成熟，是目前地震層析成像的主要方法。但是射線理論只適用于波速在一個(gè)波長(zhǎng)范圍內(nèi)變化很小的場(chǎng)合，是波動(dòng)方程的高頻近似，因此它有一定的局限性。而基于波動(dòng)方程的層析成像方法由于需要超大規(guī)模的三維數(shù)值計(jì)算，目前還有許多問題沒有解決。但波動(dòng)方程包含了地震波場(chǎng)的全部信息，比僅利用走時(shí)資料的射線追蹤層析成像更能客觀地反映地下結(jié)構(gòu)的信息，因此是未來地震層析成像的主要發(fā)展方向。層析成像技術(shù)能以圖像的方式直觀清晰地顯示地下物質(zhì)結(jié)構(gòu)的屬性，所以這種方法一產(chǎn)生就受到了極大關(guān)注，被廣泛應(yīng)用于內(nèi)部地球物理和地球動(dòng)力學(xué)、能源勘探開發(fā)、工程和災(zāi)害地質(zhì)、金屬礦勘探等領(lǐng)域。地震層析成像技術(shù)起源于20世紀(jì)30年代，自該技術(shù)應(yīng)用以來，已取得了很多重大的成果。如以美國哈佛大學(xué)和加州理工學(xué)院為代表所做的全球三維層析成像工作，首次為人類提供了地球內(nèi)部的三維結(jié)構(gòu)影像圖，其中最重要的結(jié)果是地震波速度成像結(jié)果與大地水準(zhǔn)面的相關(guān)性，地球動(dòng)力學(xué)對(duì)其給出了很好的解釋，為板塊運(yùn)動(dòng)的熱對(duì)流學(xué)說提供了證據(jù)。再如，用層析成像方法人們首次發(fā)現(xiàn)非洲超級(jí)地幔柱等大型地幔柱均起源于地幔邊界。在大洋洋脊、板塊消減帶、克拉通地區(qū)，地殼和上地幔中的火山、地殼和地幔頂部、造山帶、斷裂區(qū)和震源區(qū)等地方層析成像技術(shù)也都有大量的應(yīng)用成果。無論是能源和礦產(chǎn)等資源勘探，還是地球內(nèi)部結(jié)構(gòu)及地球動(dòng)力學(xué)研究，地震層析成像技術(shù)都是有效的、重要的技術(shù)之一。1地震層析成像技術(shù)地震層析成像涉及3個(gè)方面：數(shù)據(jù)采集、數(shù)據(jù)處理（數(shù)據(jù)正反演計(jì)算和圖像重建）、成像結(jié)果解釋。地震層析成像是采集數(shù)據(jù)的主要目的、數(shù)據(jù)解釋的基礎(chǔ)和數(shù)據(jù)處理的主要部分。地震層析成像主要包括以下幾部分：模型的參數(shù)化、正演計(jì)算地下介質(zhì)屬性的理論值（射線追蹤、波形擬和）、反演及圖像重建、反演結(jié)果的評(píng)價(jià)（分辨率分析）?，F(xiàn)分別就用于這4個(gè)方面的各種方法作一闡述。1.1建模及模型的參數(shù)化層析成像的結(jié)果是在初始模型的基礎(chǔ)上迭代反演得來的，因此初始模型與真實(shí)地下結(jié)構(gòu)接近程度直接關(guān)系到成像的結(jié)果能否準(zhǔn)確反映客觀物質(zhì)屬性。如何合理、準(zhǔn)確地描述初始模型至關(guān)重要。早期研究一般都是假設(shè)模型為均勻?qū)訝钏礁飨蛲越橘|(zhì)模型，這只是一個(gè)粗略的模型，遠(yuǎn)遠(yuǎn)不能滿足實(shí)際應(yīng)用需要。隨著研究的深入，模型逐步過渡到三維非均勻各向異性任意界面介質(zhì)模型。國際上一些標(biāo)準(zhǔn)的模型有二維的Marmousi模型，三維的鹽丘模型和逆掩模型等。Gjoystdal（1985）提出的模型生成技術(shù)，可方便地生成任意復(fù)雜結(jié)構(gòu)的殼幔模型。塊狀建模的方法，對(duì)地質(zhì)體的描述采取體f塊f面f點(diǎn)f坐標(biāo)的層次結(jié)構(gòu)。改變了傳統(tǒng)的層狀地層的建模方式，引進(jìn)了塊狀結(jié)構(gòu)的描述，并采用三角形面片構(gòu)造塊狀模型的界面，可以適應(yīng)非常復(fù)雜的三維介質(zhì)。新的建模方式從根本上改變了層狀地層建模不能適應(yīng)復(fù)雜地質(zhì)結(jié)構(gòu)的狀況。在地震層析成像技術(shù)中，由于最終反演的地下介質(zhì)屬性是通過將研究區(qū)域劃分成不重疊的多個(gè)像元，依據(jù)各像元的灰度（反演得到的地下介質(zhì)屬性）來成圖的，所以在地震層析成像中多采用網(wǎng)格的方法來進(jìn)行模型參數(shù)化。網(wǎng)格化方式也由最初網(wǎng)格內(nèi)速度均勻分布模型發(fā)展到后來的給出節(jié)點(diǎn)速度值，采用插值的方法求得網(wǎng)格內(nèi)各點(diǎn)的速度；由規(guī)則均勻網(wǎng)格發(fā)展到動(dòng)態(tài)變尺度的不規(guī)則網(wǎng)格。有人在2001年提出了一種模型參數(shù)化方法，把模型網(wǎng)格分為正演網(wǎng)格、地震網(wǎng)格（實(shí)際地下結(jié)構(gòu)模型）、反演網(wǎng)格，各種網(wǎng)格劃分密度和大小（根據(jù)射線分布情況）不同。對(duì)正演網(wǎng)格采用精細(xì)劃分，而對(duì)成像網(wǎng)格則采用相對(duì)較粗糙的劃分，這種模型參數(shù)化方法在不增加計(jì)算復(fù)雜度的情況下，提高了成像的分辨率，是一種理想的模型參數(shù)化方法。類似的研究有交錯(cuò)網(wǎng)格法，不同尺度的成像網(wǎng)格和射線追蹤網(wǎng)格，彼此通過雙曲線插值映射相互關(guān)聯(lián)。但他們采用尺寸相等的網(wǎng)格劃分方法，對(duì)于射線分布不均勻的情況有可能造成某些網(wǎng)格沒有射線通過的問題，雖然通過更小尺度的劃分可以避免上述問題，但是這無疑增加了計(jì)算量，這是和采用不等尺寸網(wǎng)格劃分方法相比不足的地方。在正演數(shù)值模擬之前，還需要做的一項(xiàng)重要工作就是數(shù)據(jù)預(yù)處理。地震層析成像結(jié)果的優(yōu)劣除了跟初始模型的選取有關(guān)，很大程度上還取決于數(shù)據(jù)空間的完備程度。如數(shù)據(jù)量的大小，數(shù)據(jù)的精度，射線分布的均勻程度及密度等。這些對(duì)于人工地震資料來說，炮點(diǎn)和接收點(diǎn)是可以人為選擇的，因此上述要求是可以得到滿足的。然而對(duì)于天然地震資料來說，只能通過數(shù)據(jù)預(yù)處理盡可能地提高成像精度，如震源深度校正、地震重新定位、時(shí)差校正、遠(yuǎn)震的高度校正和地球橢圓扁率的校正等。1.2正演數(shù)值模擬正演計(jì)算在層析成像中起著極其重要的作用。正演計(jì)算的精度和計(jì)算速度，直接決定著成像的分辨率和可靠程度。正演數(shù)字模擬技術(shù)分為求解偏微分方程的波動(dòng)方程數(shù)值模擬和由積分方程以求解波場(chǎng)傳播旅行時(shí)為主的射線追蹤數(shù)值模擬。下面分別就這兩種數(shù)值模擬技術(shù)的主要方法進(jìn)行評(píng)述。（1）射線追蹤數(shù)值模擬方法射線追蹤的方法種類較多。經(jīng)典的方法是基于初值問題的試射法和基于邊值問題的彎曲法。經(jīng)典方法存在的不足有：難以處理介質(zhì)中較強(qiáng)的速度變化，難以求出多值走時(shí)中的全局最小走時(shí)，計(jì)算效率較低。而且，試射法不能對(duì)首波和陰影區(qū)內(nèi)（射線理論不成立）的射線路徑進(jìn)行追蹤；彎曲法對(duì)于兩點(diǎn)距離較遠(yuǎn)的情況效率較低。隨著射線追蹤方法的發(fā)展，出現(xiàn)了大量不同于傳統(tǒng)方法的新型算法。這些方法的主要特點(diǎn)在于不再局限于地震波的射線路徑描述，而是直接從惠更斯原理或費(fèi)馬原理出發(fā)，采用等價(jià)的波前描述地震波場(chǎng)的特征。1）有限差分求解程函方程法Vidale基于擴(kuò)張波前的思想提出了用有限差分法求解程函方程來進(jìn)行射線追蹤的方法，開辟了一條射線追蹤的新途徑，后又于1990年將該方法推廣到三維。但是他的方法仍然沒有解決首波的射線追蹤問題，同時(shí)當(dāng)介質(zhì)中存在較大的速度間斷面時(shí)會(huì)出現(xiàn)不穩(wěn)定（不滿足因果性關(guān)系，出現(xiàn)負(fù)數(shù)開平方）。將地震波的傳播路徑近似看成一條射線的情況下，波的旅行時(shí)只與沿著波射線的速度分布有關(guān)。由于低頻波的波前恢復(fù)效應(yīng)，實(shí)際情況是當(dāng)波穿過的介質(zhì)遠(yuǎn)遠(yuǎn)小于菲涅爾帶時(shí)，這種近似的影響可以忽略，然而當(dāng)波穿過的介質(zhì)大小和菲涅爾帶相當(dāng)時(shí)，旅行時(shí)和整個(gè)菲涅爾帶內(nèi)的速度分布情況相關(guān)，這時(shí)把波的傳播路徑近似成射線對(duì)旅行時(shí)所造成的誤差是不可忽略的。2）最短路徑法最短路徑法的基礎(chǔ)是費(fèi)馬原理及圖論中的最短路徑理論，是用網(wǎng)絡(luò)節(jié)點(diǎn)之間的最小旅行時(shí)連線近似地震射線路徑的。這種方法可以同時(shí)計(jì)算出從震源到達(dá)空間所有點(diǎn)的初至走時(shí)及相應(yīng)的射線路徑，并且不受射線理論的約束，準(zhǔn)確地追蹤出陰影區(qū)內(nèi)的折射波射線路徑。波速模型的復(fù)雜性與空間維數(shù)也不會(huì)影響算法的實(shí)現(xiàn)，而且所得初至走時(shí)保證了全局最小的特性。這種算法靈活高效，實(shí)用性強(qiáng)并且克服了經(jīng)典算法的缺陷，是一種較理想的算法。但是最短路徑法的精度和速度并不比傳統(tǒng)方法強(qiáng)。當(dāng)節(jié)點(diǎn)較稀時(shí)，射線常常呈之〃字形路徑，計(jì)算出的旅行時(shí)將比實(shí)際旅行時(shí)系統(tǒng)偏大，且在波傳播方向上節(jié)點(diǎn)越少，誤差越大。另一個(gè)問題是，當(dāng)網(wǎng)絡(luò)稀疏時(shí)，特別對(duì)于速度變化平緩的區(qū)域，在兩個(gè)點(diǎn)之間常常會(huì)有幾條等時(shí)最短路徑，其中可能有一條能較好地近似真實(shí)射線路徑，而其它的卻不能，因此具有一定的不確定性。3）解析計(jì)算法解析的射線追蹤方法本質(zhì)上都是對(duì)射線方程解析求解，一般有以下幾類。基于費(fèi)馬原理的解析方法，即通過求解射線方程的極小值來求得射線的路徑。解析法適用范圍較小，因?yàn)閷?shí)際地質(zhì)構(gòu)造復(fù)雜，即速度分布比較復(fù)雜，而且解析法只能對(duì)少數(shù)特殊的速度分布實(shí)現(xiàn)射線追蹤（如速度或慢度平方是常梯度，以及慢度平方是多項(xiàng)式的情形）。4）傳統(tǒng)方法的改進(jìn)方法擬彎曲法在模型中引入了復(fù)雜形狀的速度間斷面，使用了可以計(jì)算含有速度間斷面的非均勻模型中的地震波走時(shí)和射線路徑的射線追蹤技術(shù)，適合復(fù)雜結(jié)構(gòu)地區(qū)的地震成像。該方法對(duì)射線路徑的擾動(dòng)主要分兩種情況處理。當(dāng)擾動(dòng)點(diǎn)在間斷面上時(shí)，根據(jù)斯奈爾定律，采用兩分法逐漸縮小范圍，找出在間斷面上的折射點(diǎn)：當(dāng)擾動(dòng)點(diǎn)在連續(xù)介質(zhì)內(nèi)部時(shí)，以射線路徑上3個(gè)相鄰點(diǎn)為例，先固定不相鄰兩點(diǎn)的位置，擾動(dòng)中間點(diǎn)的位置，以尋求不相鄰兩點(diǎn)之間最小走時(shí)路徑。由于這些特點(diǎn)，該算法近年來在國內(nèi)外獲得了較廣泛的應(yīng)用。5）其他方法旅行時(shí)線性插值法是Vidale差分法的一種高級(jí)形式。為了提高地震波旅行時(shí)的計(jì)算精度，提出了一種改進(jìn)的射線追蹤法一一旅行時(shí)二次/線性聯(lián)合插值法（QLTI），即在震源的近場(chǎng)采用二次插值計(jì)算旅行時(shí)，在遠(yuǎn)場(chǎng)仍用LTI算法。QLTI算法改善了近場(chǎng)旅行時(shí)的精度，降低了累加誤差，從而提高了全場(chǎng)（尤其是遠(yuǎn)場(chǎng)）旅行時(shí)的計(jì)算精度，QLTI算法較傳統(tǒng)的LTI明顯提高了精度?；诨莞乖砗唾M(fèi)馬原理求取地震波走時(shí)及其反射波射線路徑的新方法。該方法具有原理簡(jiǎn)單、易于實(shí)現(xiàn)、能適應(yīng)較為復(fù)雜地質(zhì)模型以及易于將其推廣到各向異性介質(zhì)等優(yōu)點(diǎn)，克服基本算法速度較慢的缺陷，是一種地震波走時(shí)和反射波射線路徑計(jì)算的改進(jìn)方法。在保證精度的條件下，該改進(jìn)算法的計(jì)算速度顯著提高。地震波走時(shí)層析成像方法一交錯(cuò)網(wǎng)格法，即利用高密度網(wǎng)格進(jìn)行射線追蹤，以適應(yīng)實(shí)測(cè)射線數(shù)的剖分網(wǎng)格進(jìn)行成像計(jì)算，并采用任務(wù)并行化逐炮點(diǎn)進(jìn)行射線追蹤。交錯(cuò)網(wǎng)格法的成像網(wǎng)格單元內(nèi)射線為曲線，具有較高的成像精度。方法原理簡(jiǎn)單，易于實(shí)現(xiàn)。（2）波形擬和法基于波動(dòng)方程的層析成像一般有理論地震圖法和接收函數(shù)法。由于波動(dòng)方程數(shù)值模擬實(shí)質(zhì)是求解地震波波動(dòng)方程，因此模擬的地震波場(chǎng)含了地震波的所有信息，但由于基于波動(dòng)方程的層析成像方法需要超大規(guī)模的三維數(shù)值計(jì)算，所以計(jì)算速度相對(duì)于幾何射線法要慢，且易引進(jìn)干擾波，目前還有許多困難問題沒有解決。但波動(dòng)方程包含了地震波場(chǎng)的全部信息，比僅利用走時(shí)資料僅用于模擬波的運(yùn)動(dòng)學(xué)特征的射線追蹤層析成像更能客觀地反映地下結(jié)構(gòu)的信息，因此對(duì)于研究復(fù)雜條件下的各種波場(chǎng)最為有效，具有廣闊的發(fā)展前景。目前常用的方法有偽譜法、有限元法、有限差分方法。偽譜法處理邊界靈活，是有限差分法近似階數(shù)趨于無限時(shí)的極限，它用快速傅氏變換來計(jì)算空間導(dǎo)數(shù)，計(jì)算精度要高于有限差分法。但是和有限差分方法一樣計(jì)算量大，效率較低：有限元法由于剖分的任意性及它所依據(jù)的變分原理，對(duì)含有多種介質(zhì)和自然邊界條件的處理非常方便有效，已成為解決地震波傳播數(shù)值模擬的一種重要方法。它是目前為止最精確的一種正演模擬方法，但計(jì)算量大。有限元法的主要優(yōu)點(diǎn)是適宜于模擬任意地質(zhì)體形態(tài)，可以任意三角形逼近地層界面，保證復(fù)雜地層形態(tài)模擬的逼真性。有限差分法和有限元法的主要缺點(diǎn)在于對(duì)高頻分辨的限制，對(duì)地震勘探中典型的速度和頻率，計(jì)算中需要大量的網(wǎng)格點(diǎn)，而偽譜法則相對(duì)更有效。1.3反演及圖像重建技術(shù)層析成像中的反演方法可分為線性方法和非線性方法兩種。目前非線性反演方法主要有速傳算法、模擬退火法和神經(jīng)網(wǎng)絡(luò)法等。在體波層析成像中，使用線性反演方法的較多，如奇異值分解法、共扼梯度法和最小二乘法等。最小二乘法中加上阻尼得到了阻尼最小二乘法，其實(shí)質(zhì)就是用三角矩陣分解法加上阻尼最小二乘法的超定線性方程組求解。地球物理問題大多是高度非線性問題，如層析成像問題就是一個(gè)典型的非線性問題。非線性反演方法全局搜索，且不依賴于初始模型，適用于對(duì)被研究區(qū)域的初始信息了解較少的研究，可以把其反演的結(jié)果作為初始模型進(jìn)行局部最優(yōu)圖像重建。雖然計(jì)算速度較慢，但是對(duì)于復(fù)雜的非線性反演問題效果顯著。而線性反演方法大多是人為的使非線性問題線性化，存在較大的不穩(wěn)定性，且依賴于初始模型，但是其計(jì)算速度較快。因此，在層析成像中無論是線性反演方法還是非線性反演方法，都有較多的應(yīng)用。1.4分辨率分析（反演結(jié)果的評(píng)價(jià)）層析成像解的評(píng)價(jià)是層析成像研究的重要組成部分。通過對(duì)解的分析，我們可以了解結(jié)果的可靠性、分辨率及誤差等重要信息，常用的方法有：（1）射線密度法。通過衡量每個(gè)節(jié)點(diǎn)附近的射線數(shù)量作為解的可靠性的一種評(píng)價(jià)。但這種評(píng)價(jià)方法僅給出解的可靠性的初步度量，對(duì)解的分辨率的評(píng)價(jià)還需要進(jìn)一步分析。（2）線性反演理論方法。該方法是用廣義線性反演理論，利用模型分辨率矩陣、數(shù)據(jù)分辨率矩陣和協(xié)方差矩陣來描述解的評(píng)價(jià)方法。（3）尖峰試驗(yàn)法。該方法是通過使用合成數(shù)據(jù)去獲得分辨率矩陣的列矢量，以測(cè)試方程組的病態(tài)對(duì)解的歪曲效應(yīng)。但是，它只能估計(jì)解對(duì)單個(gè)參數(shù)點(diǎn)的分辨能力，無法評(píng)價(jià)整個(gè)解的可靠性。該方法主要目的是研究已知資料的某種形狀的異常是否可以分辨。這種試驗(yàn)可以提供有關(guān)短波異常圖像的成像能力，從而有助于區(qū)分垂向分辨率和橫向分辨率的優(yōu)劣，還可以對(duì)不同大小、不同形狀異常體進(jìn)行尖峰試驗(yàn)，以便來檢驗(yàn)算法和數(shù)據(jù)對(duì)這種異常體的成像能力。（4）棋盤分辨率試驗(yàn)法。該方法的基本原理是，首先用一個(gè)人工合成數(shù)據(jù)集代替已有的觀測(cè)數(shù)據(jù)集。合成數(shù)據(jù)集由在一個(gè)特定的三維速度模型下，應(yīng)用真實(shí)的射線分布計(jì)算得到的理論走時(shí)值構(gòu)成。這一特定三維網(wǎng)格模型的速度分布（即棋盤格）是在初始一維速度模型基礎(chǔ)上，加上規(guī)則分布的擾動(dòng)值所構(gòu)成（即各節(jié)點(diǎn)的擾動(dòng)值大小相同，但正負(fù)相間地順序排列，這樣做的目的是便于分析）。然后，對(duì)合成數(shù)據(jù)集作反演計(jì)算，并把反演結(jié)果的三維速度結(jié)構(gòu)與檢測(cè)板的相似程度作為解的可靠性和分辨率的估計(jì)。由于棋盤格實(shí)驗(yàn)方法直觀、實(shí)用，便于分析對(duì)比，現(xiàn)今大部分層析成像的結(jié)果都應(yīng)用此方法進(jìn)行評(píng)價(jià)。（5）恢復(fù)分辨率實(shí)驗(yàn)法。其基本原理是，應(yīng)用反演得到的結(jié)果作為人工合成模型，在此模型中進(jìn)行射線追蹤，計(jì)算走時(shí)，同時(shí)在這個(gè)數(shù)據(jù)中加入與真實(shí)數(shù)據(jù)誤差相同數(shù)量級(jí)的隨機(jī)誤差，得到人工合成數(shù)據(jù)集。反演此數(shù)據(jù)集，得到的結(jié)果再與真實(shí)結(jié)果進(jìn)行對(duì)比，以便分析圖像恢復(fù)的情況。棋盤格實(shí)驗(yàn)設(shè)置的規(guī)則擾動(dòng)不能模擬復(fù)雜模型情況下的圖像恢復(fù)情況，相對(duì)而言，恢復(fù)分辨率實(shí)驗(yàn)要比棋盤格實(shí)驗(yàn)更接近真實(shí)情況，但不易分析和對(duì)比。2.地震層析成像技術(shù)的未來發(fā)展趨勢(shì)在以往對(duì)中國及鄰近區(qū)域的地殼上地慢三維結(jié)構(gòu)的研究中，人們比較注重地殼上地慢三維速度結(jié)構(gòu)的研究。一方面是因?yàn)橛盟俣冉Y(jié)構(gòu)可以解釋地球內(nèi)部三維結(jié)構(gòu)某些方面的特征，另一方面是因?yàn)橛?jì)算速度參數(shù)相對(duì)而言涉及的因素比較少。但是，從總的方面來看，速度參數(shù)只是利用了地震波的運(yùn)動(dòng)學(xué)信息，而忽略了地震波動(dòng)力學(xué)的特征信息。地震層析成像研究大多根據(jù)地震記錄局部上的單一觀測(cè)值反演單一物理量，方法各自獨(dú)立，表現(xiàn)出全局性和系統(tǒng)性的不足，在相當(dāng)大的程度上阻礙了地震層析成像方法的使用。地震波動(dòng)是巖石物性的綜合響應(yīng)，巖石物性參數(shù)間存在著必然聯(lián)系，單一參數(shù)反演方法對(duì)多參數(shù)多分量綜合響應(yīng)的分析顯得非常無力，因此，多參數(shù)多分量同步反演方法的研究非常必要。在以往的地震層析成像方法中，走時(shí)信息和振幅信息分別反演巖石的波速與衰減。不論走時(shí)和振幅都只能代表地震波動(dòng)的局部特征，而地震波形的整體則是地下地質(zhì)狀況和巖石物性的綜合響應(yīng)。多分量波動(dòng)層析成像突破以往只提取局部信息的作法，利用地震波動(dòng)的整體信息，通過對(duì)理論地震波動(dòng)數(shù)據(jù)與實(shí)測(cè)地震波動(dòng)數(shù)據(jù)之間波形殘差的量化分析，修改彈性波方程中的物理參數(shù)，進(jìn)而修改地質(zhì)模型的巖石物性參數(shù)。通過這樣的研究，地震層析成像研究方法將發(fā)生實(shí)質(zhì)性變化，它將使以地震波局部特征為主的研究轉(zhuǎn)變?yōu)橐缘卣鸩▓?chǎng)整體動(dòng)力學(xué)特征為主的研究。地震走時(shí)層析成像內(nèi)容：程序解釋(1)正演程序#include"stdio.h"#include"math.h"structdata{doublebegin;doubleend;doubleslope;doublelength;doubletransform[12][9];doubletransform2[108];intn;doubletime;structpoint{doublex;doubley;}point[100];}ray[144];結(jié)構(gòu)體ray[i]為第i條射線的信息：begin射線的起點(diǎn)縱坐標(biāo)end射線的終點(diǎn)縱坐標(biāo)slope射線的斜率length射線的總長(zhǎng)度transform【i】為射線在每個(gè)單元格內(nèi)大小共有兩種形式二維數(shù)組為矩陣形式，一維數(shù)組是為了方便計(jì)算而把二維數(shù)組轉(zhuǎn)化成的一位數(shù)組n計(jì)算中的中間數(shù)據(jù)，物理意義是射線所穿透的橫軸的個(gè)數(shù)，為正是向下穿，為負(fù)的是向上穿time整個(gè)射線的走時(shí)point射線與縱橫軸的交點(diǎn)(X，Y)structww{doublex;doubley;}d[100];結(jié)構(gòu)體d為了給交點(diǎn)排序的中間數(shù)據(jù)main(){FILE*fp,*fp0,*fp1,*fp2,*fp3,*fp4,*fp5,*fp6;inti,j,y,m,b,e,n,q,c,k,a;doublemm,temp1,temp2;temp1temp2為為了給點(diǎn)排序的中間數(shù)據(jù)doublex[108],xc[12][9];for(i=0;i<108;i++)x[i]=1?0/3000.0;x[20]=x[29]=1.0/5000.0;x[77]=x[78]=1?0/2000?0;/*speedmatrix*/得出慢度矩陣for(i=0;i<12;i++)for(j=0;j<9;j++)xc[i][j]=x[9*i+j];fp0=fopen(〃shijimoxing.txt〃,〃w〃)；for(i=0;i<12;i++){for(j=0;j<9;j++)fprintf(fp0,"%f”,xc[i][j]);fprintf(fp0,"\n");}/*the144rays*/144條射線的起點(diǎn)終點(diǎn)長(zhǎng)度斜率的計(jì)算for(i=0;i<12;i++){for(j=0;j<12;j++){ray[12*i+j]?begin=1?5+3?0*i;ray[12*i+j].end=1.5+3.0*j;ray[12*i+j].slope=(ray[12*i+j]?begin-ray[12*i+j].end)/45;ray[12*i+j].length=sqrt(45*45+pow(fabs(ray[12*i+j].begin-ray[12*i+j].end),2));}}/*dataout(ray)*/輸出以上數(shù)據(jù)便于驗(yàn)證fp=fopen("射線.txt","w");for(i=0;i<144;i++){fprintf(fp,〃%d%f%f%f%f\n",i+1,ray[i].begin,ray[i].end,ray[i].slope,ray[i].length);}/*dataout(speed)*/輸出慢度與所給驗(yàn)證fp1=fopen("慢度.txt”,"w”);for(i=0;i<12;i++){for(j=0;j<9;j++)fprintf(fp1,"%f”,xc[i][j]);fprintf(fp1,"\n");}/**/計(jì)算中間數(shù)據(jù)輸出fp2=fopen("中間數(shù)據(jù)?txt","w”);//for(i=0;i<12;i++)//{//for(j=0;j<12;j++){//ray[12*i+j].n=j-i;//fprintf(fp2,"%d”,ray[12*i+j].n);//}//fprintf(fp2,"\n〃)；//}///*initialize*/初始化數(shù)組for(i=0;i<144;i++)for(j=0;j<12;j++)for(y=0;y<9;y++)ray[i].transform[j][y]=0;/**/for(i=0;i<144;i++){for(j=0;j<100;j++){ray[i]?point[j]?x=0;ray[i].point[j]?y=0;}}/**/for(i=0;i<144;i++){ray[i].time=0;}fp4=fopen("y,txt”,"w”);for(i=0;i<12;i++){for(j=0;j<12;j++){/**/for(q=0;q<100;q++){d[q].x=100;d[q].y=0;}/**/計(jì)算交點(diǎn)for(m=0;m<10;m++){d[m].x=m*5.0;d[m].y=-d[m].x*ray[i*12+j].slope+ray[i*12+j].begin;}n=10;e=ray[i*12+j].n;if(e>0){for(y=1;y<=e;y++){d[n].y=ray[i*12+j].begin+3.0*y-1.5;d[n].x=(3.0*y-1.5)/fabs(ray[12*i+j].slope);n++;}}if(e<0){for(y=-1;y>=e;y--){d[n].y=ray[i*12+j]?begin+3?0*y+1?5;d[n].x=-(3.0*y+1.5)/fabs(ray[12*i+j].slope);n++;}}/**/對(duì)交點(diǎn)進(jìn)行排序b=0;while(d[b].x!=100){b++;}/**/for(c=0;c<b;c++)for(k=0;k<b-1;k++)if(d[k].x>d[k+1].x){temp1=d[k].x;d[k].x=d[k+1].x;d[k+1].x=temp1;temp2=d[k].y;d[k].y=d[k+1].y;d[k+1].y=temp2;}/**/輸出交點(diǎn)b=0;while(d[b].x!=100){ray[12*i+j].point[b].x=d[b].x;ray[12*i+j]?point[b].y=d[b],y;fprintf(fp4,"%d\t%f,%f\n",12*i+j+1,ray[12*i+j].point[b].x,ray[12*i+j].point[b].y);b++;}/**/利用兩點(diǎn)間距離計(jì)算每個(gè)射線在每個(gè)單元格上的的長(zhǎng)度for(a=1;a<b;a++){mm=sqrt(pow((ray[12*i+j].point[aT].x-ray[12*i+j].point[a].x),2)+pow((ray[12*i+j].point[aT].y-ray[12*i+j].point[a].y),2));if(ray[12*i+j].n>=0)ray[12*i+j].transform[(int)(ray[12*i+j].point[a-1].y/3)][(int)(ray[12*i+j].point[a-1].x/5)]=mm;if(ray[12*i+j].n<0)ray[12*i+j]?transform[(int)(ray[12*i+j].point[a].y/3)][(int)(ray[12*i+j].point[a-1].x/5)]=mm;}/**/計(jì)算走時(shí)for(i=0;i<144;i++)(for(j=0;j<12;j++)(for(y=0;y<9;y++){ray[i].time+=ray[i].transform[j][y]*xc[j][y];}}}輸出最后正驗(yàn)結(jié)果fp3=fopen("矩陣A?txt”,"w");for(i=0;i<144;i++){for(j=0;j<12;j++){for(y=0;y<9;y++)fprintf(fp3,"%f〃,ray[i].transform[j][y]);fprintf(fp3,〃\n");}fprintf(fp3,〃\n〃)；}fp5=fopen("矩陣A2.txt”,"w”);for(i=0;i<144;i++){for(j=0;j<12;j++){for(y=0;y<9;y++)fprintf(fp5,"%f",ray[i].transform[j][y]);}fprintf(fp5,〃\n〃)；}fp6=fopen(〃走時(shí).txt〃,〃w〃)；for(i=0;i<144;i++){fprintf(fp6,"%f\n”,ray[i].time);(2)反演程序#include"stdio,h”#include"math.h"structdata{floatA[108];floatb;floatr;}ray[144];結(jié)構(gòu)體ray【i】表示b每條射線的走時(shí)A每條射線在每條單元格內(nèi)的長(zhǎng)度r反演每條射線計(jì)算出的誤差main(){intm[108];floatx[108],R,a1,a2;x[i]反演的每個(gè)單元格的慢度intn=0,j,i,j1;FILE*fp1,*fp2,*fp3;/*,*fp3,*fp4*//**/讀取正演的結(jié)果fp1=fopen("矩陣A2?txt","r");fp2=fopen("走時(shí)?txt”,"r”);for(i=0;i<144;i++)for(j=0;j<108;j++)fscanf(fp1,〃%f〃,&ray[i].A[j]);for(i=0;i<144;i++)fscanf(fp2,"%f”,&ray[i].b);for(i=0;i<108;i++)x[i]=0;給第一次的慢度全部付零值for(i=0;i<144;i++)for(j=0;j<108;j++)printf("%f”,ray[i].A[j]);printf(〃\n〃)；for(j1=0;j1<108;j1++)算出系數(shù)矩陣中每一列不為零值的數(shù)的個(gè)數(shù){m[j1]=0;x[j1]=0;for(i=0;i<144;i++)if(ray[i].A[j1]!=0)m[j1]+=1;printf("%d”,m[j1]);n=20;迭代次數(shù)20次while(n)(誤差的計(jì)算for(i=0;i<144;i++){R=0;for(j=0;j<108;j++)R+=ray[i]?A[j]*x[j];ray[i].r=ray[i].b-R;}for(j1=0;j1<108;j1++){{a2=0;/**/利用代數(shù)重建算法迭代慢度矩陣for(i=0;i<144;i++){a1=0.0;for(j=0;j<108;j++)a1+=ray[i].A[j]*ray[i].A[j];a2+=ray[i].r*ray[i].A[j1]/a1;}x[j1]+=a2/m[j1];}//printf("%d%d\n”,j1,m[j1]);輸出中間數(shù)據(jù)}n--;}//for(j=0;j<108;j++)//{if(x[j]>2)x[j]=10;//if(x[j]<0)x[j]=0;//}輸出最后數(shù)據(jù)fp3=fopen("out.txt","w");n=0;for(i=0;i<12;i++){for(j=0;j<9;j++){fprintf(fp3,"%f”,x[n]);n++;}fprintf(fp3,〃\n〃)；}}

可編譯程序源代碼及附帶數(shù)據(jù);（3）實(shí)驗(yàn)結(jié)果，圖及文字說明（至少包含五幅圖）原始數(shù)據(jù)20次迭代200次迭代2000次迭代通過上圖正演和反演所得出的圖可以看出迭代次數(shù)越多反演得出的數(shù)據(jù)和原始數(shù)據(jù)就越接近（如右上角的異常邊緣可以看出2000次迭代的結(jié)果明顯比二十次的要好，黑色的要少）但是無論迭代多少次圖像的邊緣都會(huì)產(chǎn)生很大的異常，這應(yīng)該是