擴邊尺度對重力異常分層分離處理的影響
——以插值切割法為例

徐如剛 梁 霄 孫鴻博 儲 飛

1 安徽省地震局，合肥市長江西路558號，230031

在位場數據處理與轉換過程中，數據區域邊界往往易產生畸變等邊界效應，給資料解釋工作帶來不利影響。為盡量減少甚至消除邊界效應，保證處理后邊界數據的真實性和完整性，以便利用實測位場數據進行地質解釋，通常在處理前需要對位場數據進行必要的擴邊處理。

擴邊對位場數據處理造成的影響往往不可忽視[1-5]，且擴邊后的數據質量直接影響到處理后數據的精度。為此，很多學者圍繞如何減少邊界效應、提升處理結果可靠性等開展大量工作，并提出許多位場數據擴邊方法，如余弦擴邊方法[1]、區域場擴邊方法[1]、最小曲率擴邊方法[6]、三向擴邊方法[7]、泛克里格擴邊方法[8]、凸集投影(projection onto convex sets，POCS)擴邊方法[9]等。但這些研究主要集中于擴邊方法本身的效果與精度，而基于特定方法開展擴邊的尺度對位場數據處理影響的研究卻鮮有涉及。

為此，本文從擴邊尺度對位場數據處理影響的角度展開研究，以線性三角插值法作為位場數據的擴邊方法[10]，對傳統的四點切割算子和動態改進型插值切割算子進行改進[11]，并基于改進后的完全圓周對稱的動態改進型插值切割算子對理想模型在地表重力異常在不同深度層源進行切割分離。基于分離后的結果，從不同深度層源重力異常理論值和重力場源分離計算值之間的差異、擴邊尺度與重力場源分離誤差等角度研究擴邊尺度對重力異常分層分離結果的影響，為后續位場數據處理過程中擴邊尺度的合理選取提供參考依據。

1 插值切割分層分離原理及其技術流程

1.1 插值切割分層分離原理及改進

插值切割法是在多次切割法的基礎上發展而來的空間域劃分區域和局部場的方法[12-14]，通過引入半二階差分量，可顯著提升對異常特征的分辨能力[13]，被認為是一種較好的位場分離技術。隨著大深度位場向下延拓技術的發展，插值切割法在位場數據分離處理中得到廣泛應用[15-18]。徐世浙等[17-18]發現不同深度地層的位場效應與切割半徑之間關系密切，認為深度r以上地層的位場效應可以近似用切割半徑r切割到的局部場表示，通過r1、r2兩個不同半徑切割即可獲得r2-r1地層在地表的位場效應，插值切割分層分離技術即是基于該原理。

在大半徑切割時，現有的切割算子存在不足，如四點圓周平均切割算子會產生嚴重的窗型震蕩現象[19]，動態改進型切割算子所采用的數據點無法完全均勻分布于切割半徑所對應的圓周上[11]，因此本文綜合這2種切割算子的各自特點，對動態改進型切割算子進行改進，提出隨切割半徑增加而完全圓周對稱的動態改進型插值切割算子，并基于新的切割算子對理想模型在地表的重力異常進行切割分離。基于相同的理想模型，將本文提出的新切割算子的切割分離結果與文獻[11]的分離結果進行比較，兩者結果基本一致，表明改進后的算法具有可行性。

完全圓周對稱的動態改進型插值切割算子可確保用于計算的切割算子計算點均勻分布于切割半徑所對應圓周上，且計算點數目隨切割半徑的增加而均勻增加。根據線性三角插值法，對網格數據進行插值，重新計算切割半徑所對應圓周上的切割算子計算點處的位場值。基于插值后的切割半徑所對應圓周上的計算點處的位場值，根據式

(1)～(3)對位場進行切割分離[11]：

(1)

式中，A(x,y)表示圓心為(x,y)、切割半徑為r的圓周上N點位場的算術平均值，N隨著切割半徑的增加而增大，為切割半徑4倍，即N=4r，且均勻分布于切割半徑所對應的圓周上；G(x,y)、Gi(r)分別為圓心位于(x,y)處的位場值和切割半徑為r的圓周上點i處的位場值，圓周上i處的位場值則是通過對位場網格數據進行插值計算獲得；R(x,y)、S(x,y)分別為切割半徑r所對應的區域場和局部場；權系數g由圓心處和圓周上的位場值共同確定，反映兩者的線性程度，由式(2)和式(3)計算所得：

(2)

(3)

1.2 位場插值切割分層分離技術流程

根據文獻[10]中線性三角插值法對位場網格數據進行插值，計算1倍網格點距半徑圓周上所對應的切割算子計算點處的位場值，由計算點處的位場值根據式(1)～(3)計算1倍網格點距切割算子，基于1倍網格點距切割算子對位場數據進行切割分離，獲得的局部場可視為淺部1倍網格點距深度以上的重力效應。計算2倍網格點距半徑圓周上所對應的切割算子計算點處的位場值，由計算點處的位場值根據式(1)～(3)計算2倍網格點距切割算子，基于2倍網格點距切割算子對位場數據進行切割分離，獲得的局部場可視為淺部2倍網格點距深度以上的位場效應，該效應為第1、第2層地層的綜合重力效應，則第2層地層的重力效應為第1、第2層的綜合效應與第1層地層重力效應之差。依次類推，逐漸增加切割半徑，完成不同深度地層重力效應的分離[15,17]。

2 理想模型實驗

2.1 模型設計參數

本文以文獻[11]中給出的三維組合模型作為此次計算的理想模型，如圖1所示，模型在25 km×25 km×5 km空間范圍內布設球體、橢球體和棱柱體3種不同類型的密度異常體，模型設計參數見表1。圖2為研究區理想模型重力異常在地表的空間分布，理想模型重力異常在地表25 km×25 km范圍內以網格點距為100 m、250×250的網格狀點陣分布，模型重力異常值由Encom ModelVision軟件計算獲得。

表1 三維理想模型設計參數

圖1 理想模型空間分布Fig.1 Space distribution of ideal model

圖2 理想模型重力異常分布Fig.2 Gravity anomalies distribution of ideal model

2.2 理想模型計算

為了更直觀地認識擴邊尺度對重力異常分層分離處理的影響，本文基于改進后的完全圓周對稱的動態改進型插值切割算子，對理想模型進行實驗，以期為重力異常分層分離時擴邊尺度的合理設置提供參考依據。

2.2.1 理想模型不同深度層源的理論重力異常值計算

為確保理想模型切割分離出的不同深度層源在地表的重力異常值不產生邊界效應，在地表重力異常切割分離前，以重力異常數據格網點距的整數倍為長度對研究區重力異常數據進行擴邊。理想模型中各個異常體的埋設深度為0～5 km，模型重力異常數據的格網點距為100 m，若設定模型最大分離切割深度取50倍格網點距，則模型的最大切割深度為5 km，同時將各深度層源的切割迭代次數固定為20次，則需將研究區重力異常數據向四周擴充50×20個格網點距。如果理想模型在地表的重力異常數據擴邊尺度足夠，切割分離出的不同深度層源的重力異常則不會產生邊界效應。基于理想模型分離的最大切割深度及其在地表重力異常數據的格網點距，本文對模型重力異常數據進行擴邊，擴邊后的數據為2 250×2 250點陣，擴邊前后的數據區域范圍如圖3所示(紅色為研究區，藍色為擴邊區)。擴邊后區域范圍的地表重力異常通過Encom ModelVision軟件計算獲得，空間分布如圖4所示。分別以1～50倍格網點距為切割半徑，根據§1.2的位場插值切割分層分離技術流程對擴邊后的重力異常數據進行切割分離，并將獲得的地表重力異常稱為理想模型在不同深度層源的理論重力異常值。

圖3 研究區及其擴邊后的區域分布Fig.3 Distribution of the study area and its expanded area

圖4 擴邊后的重力異常分布Fig.4 Gravity anomalies distribution after edge expansion

2.2.2 理想模型重力異常值與重力場源分離計算值對比

由前文可知，若擴邊尺度足夠，則切割分離出的不同深度層源在地表的重力異常值不會產生邊界效應，否則，切割分離出的不同深度層源在地表的重力異常結果會包含邊界效應。本文將含有邊界效應的不同深度層源分離結果稱為理想模型重力場源分離計算值，以相同的格網點距對研究區分別進行50、100和150倍格網點距擴邊，同時基于線性三角插值法分別計算擴邊尺度為50、100和150倍格網點距的區域范圍地表重力異常。

根據分離技術流程，分別以1～50倍格網點距為切割半徑對擴邊尺度為50倍格網點距的理想模型重力異常進行切割分離，獲得1～50倍格網點距的不同深度層源在地表的重力異常。若以理想模型的重力異常值作為參考，則可獲得研究區不同深度層源在地表的重力異常理論值與重力場源分離計算值的差異，圖5為理想模型擴邊尺度為50倍格網點距的重力異常理論值與重力場源分離計算值差異的空間分布。

圖5 擴邊尺度為50倍格網點距的不同深度層源重力異常差異分布Fig.5 The distribution of gravity anomalies at different depths of layers with an expanded edge size of 50 times grid spacing

同理，計算獲得擴邊尺度為100倍和150倍格網點距的不同切割半徑對應的深度層源在地表的重力異常差異，結果如圖6和7所示。

圖6 擴邊尺度為100倍格網點距的不同深度層源重力異常差異分布Fig.6 The distribution of gravity anomalies at different depths of layers with an expanded edge size of 100 times grid spacing

圖7 擴邊尺度為150倍格網點距的不同深度層源重力異常差異分布Fig.7 The distribution of gravity anomalies at different depths of layers with an expanded edge size of 150 times grid spacing

由圖5～7可以發現：1)理想模型重力異常分層分離結果中的邊界效應與模型研究區數據的擴邊尺度具有密切關系，若模型分層分離的深度相同，隨著模型數據擴邊尺度的增加，分層分離結果受到邊界效應的影響則越弱；2)理想模型重力異常分離結果中的邊界效應與其切割分離的深度具有密切關系，若模型研究區擴邊尺度相同，切割分離的深度越深，分離后的結果中包含的邊界效應則越明顯；3)邊界效應與異常體的規模也有一定關系，若模型研究區數據的擴邊尺度、分層分離深度相同，異常體的規模越大，分層分離結果中的邊界效應則越明顯，擴邊尺度需進一步增加，以減弱分層分離結果中的邊界效應；4)若模型數據擴邊尺度固定，分層分離結果中邊界效應的空間展布不是由數據區域邊界向中心區域均勻分布，而與異常體的規模和空間位置有關，異常的規模和幅度越大，邊界效應則越明顯。

2.2.3 擴邊尺度與重力場源分離誤差分析

根據均方誤差公式對擴邊尺度與重力場源分離誤差進行分析，分別統計理想模型在50、100和150倍格網點距擴邊尺度條件下1～50倍格網點距對應的不同深度層源在地表的重力異常理論值與重力場源分離計算值差異的均方誤差：

(4)

式中，M、P分別為理想模型重力異常格網的線數及其點數，vt、vc分別為理想模型不同深度層源在地表的重力異常理論值和重力場源分離計算值。δ計算結果如圖8所示。

圖8 不同擴邊尺度、不同深度層源重力異常差異均方誤差Fig.8 The mean square error of gravity anomalies difference at different expanded edge sizes and depths

圖8中橫坐標為理想模型切割所選取的切割半徑，最大切割半徑為50倍格網點距；縱坐標為不同切割半徑所對應的重力異常差異的均方誤差。從圖8可以看出：1)當理想模型數據的擴邊尺度一定時，隨著切割半徑的增加，即切割深度增加，不同深度層源的地表重力異常差異的均方誤差不斷增大，即分層分離結果中的邊界效應越顯著；2)當切割深度一定時，隨著數據擴邊尺度的增加，不同深度層源的重力異常差異的均方誤差迅速衰減，表明擴邊尺度越大，分層分離結果中的邊界效應顯著減弱；3)隨著理想模型數據擴邊尺度的增加，模型重力異常差異的均方誤差衰減變緩，表明擴邊尺度增加到一定程度再繼續增加時，對邊界效應的改善有限；4)理想模型數據的擴邊尺度與切割的最大深度存在近似2∶1的關系，若要使不同深度層源的重力異常分離結果的邊界效應不明顯，則理想模型數據的擴邊尺度不應小于最大切割深度的2倍；5)隨著模型數據擴邊尺度的增加，若用小尺度切割半徑進行切割時，累積誤差效應顯著增加，表明模型數據的擴邊尺度并非越大越好。

同理，分別計算理想模型數據的擴邊尺度基數為40倍格網點距(40、80、120和800)、30倍格網點距(30、60、90和600)及20倍格網點距(20、40、60和400)條件下理想模型的地表重力異常理論值與重力場源分離計算值差異的均方誤差，結果如圖9所示。

通過對比圖8和9可以發現，隨著模型數據擴邊尺度的增加，在減弱分離結果邊界效應的同時也會使分離結果中引入計算誤差。當模型數據的擴邊尺度過大，且切割半徑較小時，分離結果中計算誤差顯著，且隨著擴邊尺度增加，誤差累積效應越明顯，即數據的擴邊尺度并非越大越好，適當即可，如模型數據的擴邊尺度基于50倍和40倍格網點距(圖8和9(a))。當模型數據的擴邊尺度較小、切割分離半徑較大時，分離結果中包含的誤差累積效應顯著減弱，如模型數據擴邊尺度基于30倍和20倍格網點距(圖9(b)和9(c))。

圖9 不同擴邊尺度、不同深度層源重力異常差異的均方誤差Fig.9 The mean square error of gravity anomalies difference at different expanded edge sizes and depths

3 實際資料計算

圖10為五河周邊區域1∶20萬比例尺的布格重力異常數據，插值后的數據格網間距為1 000 m×1 000 m，本文基于改進后的插值切割算子對研究區實測布格重力異常數據進行不同深度層源的切割分離，切割分離最大深度為10倍格網間距，即10 km，計算獲得1 000 m、2 000 m、3 000 m和6 000 m不同深度層源在地表的重力異常，分離后的各深度層源在地表的重力異常如圖11所示。基于改進后的插值切割算子的分離結果與文獻[11]基本一致，且隨著切割深度的增加，兩者一致性更好，但淺部分離結果(圖11(a))與文獻[11]在異常幅度方面存在較為明顯的差異，這種差異可能是由于本文與文獻[11]使用不同分辨率的重力異常數據進行切割分離所致。

圖10 布格重力異常Fig.10 Bouguer gravity anomalies

圖11 擴邊尺度為40倍格網點距的不同深度層源重力異常在地表的空間分布Fig.11 The distribution of gravity anomalies at different depths of layers with an expanded edge size of 40 times grid spacing

根據均方誤差公式對實測布格重力異常數據的擴邊尺度與重力場源分離誤差進行分析，以40倍格網點距擴邊尺度為參考，計算10、20和30倍格網點距擴邊尺度條件下1～20倍格網點距對應的不同深度層源在地表的重力異常差異的均方誤差，結果如圖12所示。圖12與圖8總體上表現出相似的特征，同時圖12的分層分離結果更好地說明在邊緣區域存在較大異常時，需要在擴邊尺度與最大切割深度近似2∶1的基礎上進一步增加擴邊尺度，以減弱分層分離結果中的邊界效應，如擴邊尺度與最大切割深度關系近似為3∶1。

圖12 不同擴邊尺度、不同深度層源重力異常差異的均方誤差Fig.12 The mean square error of gravity anomalies difference at different expanded edge sizes and depths

4 結 語

本文基于改進后的插值切割算子，對三維組合理想模型在地表的重力異常進行不同深度層源的切割分離，從不同深度層源的重力異常差異、擴邊尺度與重力場源誤差分析方面，研究擴邊尺度對重力異常分層分離處理的影響，為位場數據處理過程中擴邊尺度的合理選取提供參考依據。

1)實驗結果表明，數據區域的邊界效應與擴邊尺度、分離切割深度、異常的規模和幅度及空間位置等因素相關。當切割分離深度、異常體的規模和異常幅度、空間位置一定時，數據擴邊尺度越大，分層分離結果中的邊界效應越弱；擴邊尺度越小，邊界效應越強。

2)擴邊尺度的選擇應綜合考慮分析資料所需滿足的精度和數據計算效率，擴邊尺度不宜過大，若擴邊尺度過大，且切割半徑較小，分離結果中計算誤差顯著，擴邊尺度與重力異常分層分離的最大切割深度至少應滿足2∶1的要求，才能明顯改善邊界效應和計算誤差。對含有較大重力異常的邊界進行切割分離時，應在2倍擴邊尺度的基礎上再進一步擴邊，基于此關系，根據所要分離的最大切割深度確定數據區域的擴邊尺度，或根據數據區域的擴邊尺度確定最大切割深度。

3)數據區域邊界效應的空間展布不是由數據區域邊界向中心區域均勻分布，而與異常體的規模和空間位置有關。

致謝:本文得到中國地震局地震研究所胡敏章博士、合肥工業大學葛粲博士、中國地震局第一監測中心劉金釗博士的支持和幫助，在此表示誠摯感謝。