λ－f域拋物Radon變換多次波壓制方法

王維紅，張 振，石 穎，李 瑩

（東北石油大學地球科學學院，黑龍江大慶 163318）

λ－f域拋物Radon變換多次波壓制方法

王維紅，張 振，石 穎，李 瑩

（東北石油大學地球科學學院，黑龍江大慶 163318）

空間截斷效應降低傳統Radon變換壓制多次波的精度.提出一種λ－f域拋物Radon變換方法，與傳統最小二乘拋物Radon變換相比，地震數據變換到λ－f域后表現為直線的映射，更便于設計濾波函數，從而有效消除傳統多次波壓制過程中空間截斷效應的影響.另外，新變量λ的引入，使Radon正反算子成為與頻率無關的矩陣，可大幅度提高計算效率.理論模型和實際海洋地震數據多次波壓制試算結果表明，提出的λ－f域拋物Radon變換法可有效實現多次波壓制，提高地震數據的信噪比，具有較強的實用性.

Radon變換；多次波壓制；λ－f域；計算效率

0 引言

多次波壓制是海洋地震資料處理中的重要預處理環節［1－4］，因為多次波的存在降低地震數據的分辨率，影響地震數據的成像質量［5－7］和后續地震解釋的精度［8］，有效地壓制多次波在地震數據處理中顯得尤為重要.Radon變換是疊前地震數據處理中壓制多次波的有效方法，人們研究Radon變換壓制多次波的方法，Hampson D［9］提出頻率域最小二乘拋物Radon變換方法，但是在Radon域存在剪刀狀“拖尾”現象，致使濾波函數設計困難，從而導致多次波壓制效果不理想.Sacchi M等［10］改進Hampson方法，提出頻率域高分辨率拋物Radon變換方法；該方法在處理簡單數據時可以有效地解決變換域“拖尾”現象，但對復雜介質的地震數據，無法得到高精度多次波壓制效果.也有學者提出時間域Radon變換方法［11－12］彌補頻率域Radon分辨率較低的缺陷，但時變矩陣的求解無法應用Levinson遞推快速計算方法，導致計算效率大幅度降低.

筆者提出一種基于λ－f域拋物Radon變換壓制多次波的方法，通過變量λ的引入，該方法的Radon變換算子與頻率無關，從而比傳統的Radon變換具有更高的計算效率［13］.同時，經λ－f域拋物Radon變換，拋物線型的地震同相軸在λ－f域變為具有一定斜率的直線，在一定程度上避免多次波壓制過程中“拖尾”現象的干擾［14］.經過理論模型及實際地震資料的地震道集數據多次波壓制試算，該方法具有較好的多次波壓制效果.

1 最小二乘拋物Radon變換

拋物Radon正反變換時間域離散積分形式為

式（1－2）中：d（x，t）為時間—炮檢距域地震數據；u（q，τ）為對應的Radon域的模型空間；q為曲率參數；τ為Radon域零偏移距截距時間；t為雙程旅行時；Nx為地震數據的道數；Nq為曲率參數q的個數；xn為偏移距.

由于時間t與截距時間τ存在線性關系，可經由頻率域高效計算拋物Radon變換，所以方程式（1）、（2）沿時間方向做傅里葉變換后，頻率域形式表示為

式（3－4）中：u（q，f）、d（x，f）分別為u（q，τ）和d（x，t）對應的傅里葉變換；f為頻率；j為虛數單位.對于單個頻率成分，方程式（4）的矩陣形式為

其中，Radon算子為

由于L與LH不是互逆的，一般不直接采用算子LH計算Radon變換，而是采用最小二乘方法定義拋物Radon正變換，即

式中：I為單位矩陣；μ通常取LHL矩陣主對角線數據的0.01倍［15］，增強矩陣求解的穩定性.

式（5）、（7）即為正反Radon變換的關系式，由于式（5）、（7）中的Radon算子與頻率有關，所以在式（5）、（7）求解過程中，每計算一個頻率分量，需要計算一次Radon算子及算子的逆，降低計算效率.

2 λ－f域拋物Radon變換

定義新變量λ＝fq，式（4）在λ－f域的Radon變換定義為

式中：M為變量λ的個數.

與式（7）對應的λ－f域拋物Radon正變換為

其中，Radon算子L（λ）表達式為

λ－f域拋物Radon變換的算子L（λ）是關于偏移距參數xn及λ的復數矩陣，與頻率無關，在整個λf域拋物Radon正反變換過程中僅需計算一次，相比于傳統τ－q域拋物Radon變換，計算效率有很大程度的提高.并且LHL為主對角線元素相等的Toeplitz矩陣，所以Levinson遞推法快速求解依然適用.

理論模擬炮記錄的Radon變換域見圖1.圖1（a）為含有3個同相軸的單炮地震記錄，共100道，道間距為20m，時間采樣間隔為4ms，每道600個采樣點.對圖1（a）的地震數據采用最小二乘τ－q域拋物Radon正變換，理論上在Radon域得到3個點，然而實際上是3個含有很大假象的剪刀狀的結果（“拖尾”現象）（見圖1（b））.這是因為在Radon變換離散運算中，地震數據的有限孔徑引起截斷效應.在Radon變換壓制多次波過程中，由于存在截斷效應的干擾，有效波與多次波在Radon域的假象可能交織在一起，使得有效波與多次波的假象得不到完全分離，進而導致波場分離不徹底.

為消除空間截斷效應的影響，對圖1（a）的地震數據進行λ－f域拋物Radon變換，得到λ－f域剖面（見圖1（c）），地震數據在經過λ－f域拋物Radon變換后，能量分布在以q為斜率的直線上，并未因空間截斷效應在λ－f域產生假象，在一定程度上避免“拖尾”現象的干擾.

3 離散采樣及濾波函數

3.1 離散采樣

為避免在λ－f域拋物Radon變換過程中出現假頻，對采樣間隔參數Δλ和λ的范圍進行合理離散采樣，其中：

式（11－12）中：λmax、λmin分別為λ的最大和最小值；xmax、xmin分別為數據的最大和最小偏移距；ΔTmax、ΔTmin分別為地震數據在時空域的最大與最小動校正時差；fmax為數據最大有效頻率.

3.2 濾波函數

根據λ定義，地震數據d在經過λ－f域拋物Radon變換后，一次波與多次波呈現過原點以不同曲率參數q為斜率的直線，引入濾波矩陣G（λ，f），對有效波和多次波進行分離.

式中：ΔT0為有效波與多次波分離所對應的動校正時差.有效波在λ－f域可表示為

變換到x－f域的有效波可表示為

4 數據試算

4.1 理論模型

為驗證λ－f域拋物Radon變換在壓制多次波中的有效性，合成含表面多次波的模型數據，見圖2（a）.理論模型數據共100道，道間距為20m，時間采樣間隔為4ms，每道601個采樣點.分段部分動校正后的剖面見圖2（b），一次波同相軸向上彎曲，曲率參數q為負值，為校正過量；多次波同相軸向下彎曲，曲率參數q為正值，為校正不足.地震數據在經過部分動校正后，有效波與多次波在變換域內得到更加有效的區分.對圖2（b）的地震數據應用文中方法得到λ－f域剖面（見圖2（c））.根據λ定義，一次波與多次波在λ－f域剖面上為經過原點且分布在原點兩側不同斜率的直線上.根據式（13）選取合適的濾波算子對圖2（c）進行多次波切除，得到濾波后λ－f域剖面（見圖2（d））.λ－f域剖面變換到時空域得到的單炮記錄見圖2（e），壓制多次波后炮記錄只有5個有效波同相軸，在圖2（a）的其他多次波得到有效壓制，即地震資料的信噪比得到有效提高.對理論模型數據應用最小二乘拋物Radon變換方法壓制多次波后的剖面見圖2（f），模擬的地震數據在1.8s處，有效波同相軸與多次波同相軸在時間域難以分開，在Radon變換域由于受空間截斷效應的影響，兩者的像在變換域也難以有效分開，所以影響濾波函數的設計，使最小二乘拋物Radon變換多次波壓制的效果不理想，即多次波有較大殘余能量，同時有效波的振幅也受到一定程度的損失.

圖2 合成炮記錄的多次波壓制效果Fig.2 Multiple suppression result on synthetic shot record

4.2 實際地震數據多次波壓制

為進一步測試文中方法在壓制多次波中的穩定性和效果，對海洋地震數據進行多次波壓制試算（見圖3）.圖3（a）為某工區海洋地震資料的一個CMP道集，數據共60道，道間距為50m，時間采樣間隔為4 ms，時間方向8×103個采樣點.圖3（a）箭頭在0.8s處的一階表面多次波較為發育，應用λ－f域拋物Radon變換方法對該地震數據進行濾波處理，得到多次波壓制后的道集記錄見圖3（b），能量很強的一階表面多次波得到有效的壓制，一次波的同相軸清晰可見，地震數據的信噪比得到有效提高，可為后續地震資料成像處理提供數據基礎，實際數據的計算表明文中方法在壓制海洋地震資料的多次波方面具有有效性和實用性.

圖3 實際數據多次波壓制結果Fig.3 Multiple suppression result on field data

5 結論

（1）λ－f域拋物Radon變換通過引入新變量λ，Radon算子的特點和傳統最小二乘方法不同，其正反算子僅需計算一次，大幅提高該方法的計算效率.

（2）拋物線型同相軸在λ－f域拋物Radon變換域為過原點的具有一定斜率的直線，在一定程度上避免傳統拋物Radon變換壓制多次波過程中空間截斷效應的影響，并且其多次波壓制濾波函數的設計也更為簡單.

（3）應用λ－f域拋物Radon變換可有效提高地震數據的信噪比，具有穩定、有效和實用的特點.

（References）：

［1］ 薛東川.幾種疊前逆時偏移成像條件的比較［J］.石油地球物理勘探，2013，48（2）：222－227.

Xue Dongchuan.Comparison of several imaging condition for prestack reverse－time migration［J］.Oil Geophysical Prospecting，2013，48（2）：222－227.

［2］ 馬繼濤，Sen M K，陳小宏.平面波域反數據處理壓制多次波方法研究［J］.地球物理學報，2009，52（3）：808－816.

Ma Jitao，Sen M K，Chen Xiaohong.Multiple attenuation using inverse data processing in the plane－wave domain［J］.Chinese Journal of Geophysics，2009，52（3）：808－816.

［3］ 劉伊克，Sun Hongchuan，常旭.基于波射線路徑偏移壓制多次波［J］.地球物理學報，2004，47（4）：697－701.

Liu Yike，Sun Hongchuan，Chang Xu.Multiple removal by wavepath migration［J］.Chinese Journal of Geophysics，2004，47（4）：697－701.

［4］ 熊登，趙偉，張劍鋒.混合域高分辨率拋物Radon變換及在衰減多次波中的應用［J］.地球物理學報，2009，52（4）：1068－1077.

Xiong Deng，Zhao Wei，Zhang Jianfeng.Hybrid－domain high－resolution parabolic Radon transform and its application to demultiple［J］.Chinese Journal of Geophysics，2009，52（4）：1068－1077.

［5］ 郭雪豹，王建民，王維紅，等.基于GPU并行加速的VSP數據逆時偏移［J］.東北石油大學學報，2014，38（2）：58－62.

Guo Xuebao，Wang Jianmin，Wang Weihong，et al.VSP data reverse－time migration based on GPU parallel acceleration［J］.Journal of Northeast Petroleum University，2014，38（2）：58－62.

［6］ 石穎，陸加敏，柯璇，等.基于GPU并行加速的疊前逆時偏移方法研究［J］.東北石油大學學報，2012，36（4）：111－115.

Shi Ying，Lu Jiamin，Ke Xuan，et al.Prestack reverse time migration based on GPU parallel accelerating algorithm［J］.Journal of Northeast Petroleum University，2012，36（4）：111－115.

［7］ 田東升，王云專，李義鵬，等.單程和雙程波動方程疊前深度偏移方法［J］.東北石油大學學報，2014，38（4）：39－44.

Tian Dongsheng，Wang Yunzhuan，Li Yipeng，et al.One－way and two－way wave equation pre－stack depth migration approaches［J］.Journal of Northeast Petroleum University，2014，38（4）：39－44.

［8］ 張進鐸.地震解釋技術現狀及發展趨勢［J］.地球物理學進展，2006，49（2）：578－587.

Zhang Jinduo.Present status and future trend of seismic data interpretation techniques［J］.Progress in Geophysics，2006，49（2）：578－587.

［9］ Hampson D.Inverse velocity stacking for multiple elimination［J］.Journal of the Canadian Society of Exploration Geophysicists，1986，22（1）：44－55.

［10］ Sacchi M，Ulrych T.High－resolution velocity gather and offset space reconstruction［J］.Geophysics，1995，60（4）：1169－1177.

［11］ Trad D，Ulrych J，Sacchi M.Accurate interpolation with high－resolution time－variant Radon transform［J］.Geophysics，2002，67（2）：644－656.

［12］ Ng M，Perz M.High resolution Radon transform in the t－x domain using"intelligent"prioritization of the Gauss－Seidel estimation sequence［J］.G xandd Abra，2004，18（3）：260－263.

［13］ 石穎，張振，李婷婷，等.λ－f域加權拋物Radon變換地震數據重建方法研究［J］.地球物理學進展，2014，29（4）：1752－1757.

Shi Ying，Zhang Zhen，Li Tingting，et al.Investigation on seismic data reconstruction by domain weighted parabolic Radon transform method［J］.Progress in Geophysics，2014，29（4）：1752－1757.

［14］ 謝俊法，孫成禹，韓文功.迭代拋物Radon變換法分離一次波與多次波［J］.石油地球物理勘探，2014，49（1）：76－81.

Xie Junfa，Sun Chengyu，Han Wengong.Iterative parabolic Radon transform for primary and multiple separations［J］.Oil Geophysical Prospecting，2014，49（1）：76－81.

［15］ Schonewille M A，Duijndam A J W.Parabolic Radon transform，sampling and efficiency［J］.Geophysics，2001，66（2）：667－678.

DOI 10.3969／j.issn.2095－4107.2015.01.003

TE132.1，P631

A

2095－4107（2015）01－0017－06

2014－12－01；編輯：陸雅玲

國家自然科學基金項目（41474118）；國家高技術研究發展計劃（863計劃）項目（2012AA061202）

王維紅（1975－），男，博士，副教授，主要從事地震資料數字處理方面的研究.