基于泊松碟采樣的小生境遺傳算法剩余靜校正

周 強 曹琳昱 陳愛萍

(①中國石油集團東方地球物理公司物探技術研究中心成都分中心，四川成都 610213； ②中國石油集團東方地球物理公司研究院西南分院，四川成都 610213； ③中國石油集團山地地震技術試驗基地，四川成都 610213)

1 引言

山地復雜構造是油氣資源的富集區之一。在其形成過程中，經歷了多期構造運動，使地表形態和地下構造復雜，靜校正問題嚴重。無論采用折射或層析方法構建的近地表結構，都只是對速度、厚度的一個近似解，不可避免存在剩余靜校正問題[1-4]。這將直接影響地震資料成像效果，從而影響對地下構造認識的精度和可信度。

目前主要的反射剩余靜校正方法有三種。第一種是高斯賽德爾分解方法，以能量最大確定延遲時，建立反射時間方程組，采用高斯賽德爾迭代求取剩余靜校正量，計算速度快。現今主流商業軟件大多配備此類模塊，該類方法適用于絕大多數地震資料。但由于采用線性方法求解剩余靜校正這種典型的非線性問題[5，6]，在處理山地低信噪比資料時，多輪之后會出現剩余靜校正量不收斂的現象。第二種是能量迭代法，它對低信噪比資料有一定適應性，但因同一物理點對應的靜校正量須滿足隨機性假設，否則存在累計誤差，進而影響計算結果。第三種是非線性方法，包括模擬退火、遺傳算法等。常采用不同搜索策略開展全局尋優[7-11]，但存在一個通病就是運算速度極慢，難以滿足實際數據處理中的高效要求。因此，亟待研發一種能全局尋優、快速解決山地低信噪比資料剩余靜校正問題的、可工業化實施的非線性剩余靜校正方法。

本文基于常規遺傳算法，利用泊松碟采樣生成初始解； 在進化過程中通過多線程構建多種群小生境，采用計算多核并行搜索最優剩余靜校正量； 最后應用理論模型和實際資料從效果和效率兩方面驗證了該方法的有效性。

2 方法技術

Wilson等[12，13]率先將遺傳算法應用于求解地震信號剩余靜校正量。以疊加能量最大為準則，對多個剩余靜校正量組合進行交換、選擇、變異操作，最終優選能量最大的靜校正量組合，有效地避免了線性方法易陷入局部極值的問題。

常規遺傳算法剩余靜校正在編代碼實現時，因C語言自帶的隨機函數生成初始種群的隨機性不夠，采用單一種群容易產生“早熟”現象； 另外，計算效率低也是其明顯缺陷。針對常規遺傳算法進行剩余靜校正存在的弊端，本文設計出一種基于泊松碟采樣的小生境遺傳算法(圖1)，用于更準確地求取地震資料剩余靜校正量。

圖1 本文方法概要

2.1 泊松碟采樣均勻初始解

隨機數通常以系統時間為種子生成，隨機剩余靜校正量基本在同一時間完成，使隨機性不夠，導致初始解空間不均勻，包含最優解的可能性低，對后期的遺傳進化依賴性強，進而影響遺傳算法的全局尋優，且搜索耗時長[14]。針對該問題，擬采用泊松碟采樣思路。

泊松碟均勻初始解空間的方法是每個剩余靜校正量隨機生成，任意兩個值的差值大于閾值，以該差值為半徑的圓內不再生成其他值。具體做法如下。

統計疊前地震數據體中炮點總數(Ns)及檢波點總數(Nr)，結合用戶設定的靜校正量范圍[a，b]，并采用實數編碼，獲取初始解空間

ri,j=U[a,b]

i∈[1,N]j∈[1,Ns+Nr]

(1)

式中：ri，j為初始解，且i為個體序號，j為基因位置；N為一個種群的個體數；U為解空間；a、b分別為最小、最大剩余靜校正量。

為了增強初始解遍布解空間的均勻度，采用泊松碟采樣方法在全部解空間中抽取相應的個體。具體方法為：首先計算解空間的規模，將解空間的最小范數解與最大范數解之間的距離作為解空間的規模參數

(2)

式中γ、κ分別為解空間中的最小范數和最大范數解。采用下列閾值

(3)

若該閾值不滿足，則重復計算；通過篩選迭代，最終獲得下列初始種群矩陣

(4)

2.2 種群內遺傳進化

種群內采用“精英+賭輪”選擇、余下個體與精英兩兩交叉、小擾動變異策略進行遺傳進化，實現全局尋優。

計算種群中每個個體對應的目標函數值，即為種群中個體的適應度值。以剖面能量為例，目標函數計算公式為

]2

si,ns,ri,nr∈R

(5)

式中：Ei為疊加剖面能量；s、r分別為炮點和檢波點的剩余靜校正量；d為疊前地震信號幅值；nr、ns分別為第i個個體對應的炮點、檢波點號；i為個體序號；y為CDP號；h為炮檢距；t為時間。根據個體適應度值，保留k個適應度最高的個體作為精英，對余下個體計算選擇概率

(6)

式中：E為個體對應的目標函數；P為選擇概率；R為個體。之后選擇直接進入子代的個體，生成隨機數i，若滿足

(7)

則選擇第i個父代個體進入子代，其余個體依據

r∈U[0,1]

(8)

實現交換策略。若

r

(9)

則有

(10)

否則不交換；若交換，即由種群中最優個體與其交換產生新個體。采取兩點交換策略，即隨機產生兩個交換位置，將兩個待交換個體中兩個交換位置之間的基因段進行交換，產生子代新個體； 否則，父代個體直接進入子代成為子代個體。

對交換后種群，保留k個適應度最高的個體作為精英，其余個體依據變異概率Pm判斷是否進行變異，即產生[0，1]間隨機數r。若r

(11)

式中：k為常數，可據實際設定；E0為初始目標函數；Eiter為當前代最優解的目標函數； 靜校正量范圍[a，b]滿足

[ri,j-c,ri,j+c]∩[a,b]

(12)

2.3 種群間聯姻

常規遺傳算法具有極大概率出現早熟現象，即種群中某單個優秀個體的基因會漸漸“統治”其他基因，造成整個種群無法跳出該個體基因的“統治”而陷入局部極值，進而導致算法無法得到最優解。

為了增加基因段的多樣性，同時減少早熟現象的發生，采取多種群同時進行遺傳操作。如圖2所示，在上一輪遺傳進化之后，進行聯姻操作。具體做法是將每個種群中的適應度函數值最大的個體的大部分基因交換至另一種群適應度函數值中等的個體，再開始下一輪自由進化。該方式下若干不同特征的種群獨立進行遺傳操作，獲得所有種群的一個最優解后，在各種群間進行通婚操作，使得全部種群都可分享優良基因。

圖2 多種群小生境遺傳算法操作示意圖

為進一步提高計算效率，調用omp_set_num_ threads函數設置進程個數。該進程個數取決于設置的小生境種群個數，通常小于CPU的核數，但一般要求不小于4，以確保種群良性發展。每個線程處理一個種群的遺傳操作，包括每個個體的基因選擇、交換機變異。采用多線程并行方式處理小生境遺傳，能充分利用計算機多核更快地尋優。

如此迭代計算，直至達到預設的迭代次數或誤差門檻值。幾個種群中疊加能量大的個體，即為估算的炮點、檢波點剩余靜校正量。

3 理論模型測試

構建一個“雙平層+單斜層”理論模型。其二維理論數據采用30Hz雷克子波合成，1ms采樣，采樣長度為1000ms； 觀測系統為1785-15-30-15-1785，120道接收，檢波點、炮點樁號范圍分別是1～318、61～259，滿覆蓋次數為30。

圖3a和圖3c分別是炮點、檢波點靜校正量：黃線為期望值，范圍是-20～20ms，有較強隨機性； 藍線是常規剩余靜校正量； 粉線是本文方法計算的小生境遺傳算法剩余靜校正量。三種顏色曲線整體上吻合較好，說明大部分剩余靜校正問題已得到解決。但從炮點、檢波點靜校正量誤差值(圖3b、圖3d)可明顯看出，藍線有個別異常值，說明常規剩余靜校正方法存在局部解，而本文提出的小生境遺傳算法剩余靜校正能全局尋優，未陷入局部解。

從疊加剖面來看，常規剩余靜校正的局部解導致剖面上同相軸出現錯斷及虛假同相軸(圖4a)現象，而本文剩余靜校正方法除了小號附加段因覆蓋參數不足出現稍許彎曲外，其余剖面都恢復了構造形態，且能量得到很好聚焦(圖4b)。

進一步分析共中心點道集(CMP1281)。將加入了隨機擾動量的道集(圖5a)、常規剩余靜校正后道集(圖5b)與用本文方法的剩余靜校正后道集(圖5c)進行對比，可見常規剩余靜校正的兩個局部解造成兩道地震道發生了較大時移，在疊加剖面上表現為同相軸的錯斷。

圖6a展示的是小生境遺傳算法和常規遺傳算法計算的疊加能量(縱軸)隨大循環(即小生境算法中的聯姻)次數收斂的情況，可見小生境遺傳算法最終的疊加能量更大，而常規算法收斂于某一能量值就不再上升，即所謂的“早熟”。

圖6b展示的是在第一次循環時，即還沒發生聯姻前，四個種群各自選擇、交叉、變異。其進化過程各不相同，即疊加能量的具體變化不同，而最終的疊加能量是A種群最大，則A 種群中能量最大的精英基因將替代其他三個種群能量處于中間值的個體的部分基因； 之后，四個種群進入下一次大的進化循環。

圖3 剩余靜校正量及其誤差

圖4 理論模型成像效果對比

圖5 共中心點道集

圖6 能量收斂圖

4 實際應用及效果

新疆某山地地區斷裂發育，地層傾角較大，地下構造復雜；近地表結構復雜，縱橫向速度、厚度變化劇烈。無論采用初至折射還是層析折射靜校正方法，都只能取得對近地表的近似成像效果，不可避免地存在較嚴重的剩余靜校正問題。

基于此，在進行多輪常規剩余靜校正和速度迭代更新后，成像品質仍不能很好地滿足地質任務要求。應用本文提出的基于泊松碟采樣的小生境遺傳算法剩余靜校正方法處理后，成像效果(同相軸連續性、波組特征)有了明顯改善。

圖7a是三維地震M測線經多輪常規剩余靜校正和速度迭代更新后的疊加剖面，存在的問題是構造主體部位同相軸連續性不好，地震反射特征雜亂，波組特征不夠清晰。應用本文方法后，剖面整體成像效果(圖7b)得到改善：同相軸變得連續；尤其是通過細節對比，構造頂部波組特征變得清晰，兩翼與頂部的銜接自然連續(圖8)。

水平疊加時，剩余靜校正和速度更新都需反復迭代，耗費大量時間。而本文方法實現了多線程編程，應用計算機多核優勢縮短了估算剩余靜校正量的時間。通過效率對比(表1)發現，達到能量增加率終止條件，常規遺傳算法需608min，而本文方法只需49min，其處理效率得到大幅度提高。

表1 效率對比表

圖7 應用本方法前(a)、后(b)三維地震M測線成像效果對比

圖8 應用本文方法前(a)、后(b)疊加剖面構造頂部細節對比

5 結論與討論

本文提出的估算剩余靜校正量的方法，目標函數為反射疊加能量，求解對象為炮點、檢波點剩余靜校正量； 在此基礎上，應用泊松碟采樣對求解的靜校正量實數編碼，增加初始解遍布解空間的均勻度； 構建多種群小生境，種群內采用“精英+賭輪”選擇、余下個體與精英兩兩交叉、小擾動變異策略進行遺傳進化，種群間采用多線程同時進化，選出幾個種群中適應度最大個體與其他種群聯姻，增加基因片段的多樣性，避免“早熟”，提高尋優的效率與精度。理論模型和實際數據測試表明，該方法可顯著且快速改善剖面成像效果，為復雜地區(尤其是低信噪比)資料提供了一種解決高頻靜校正問題的實用方法。

本文方法采用了多線程方式，效率提高了約10倍。但隨著對地震資料品質要求的提高，高覆蓋次數及新老資料的連片處理使得地震數據量急劇增加，要使本文方法具有更強適應性，還需基于地表一致性原則構建MPI框架，以增強超大數據處理能力； 同時應用多節點、多線程對剩余靜校正量進行估算，進一步提高計算速度及處理效率。