劉福香， 王萬(wàn)銀， 紀(jì)曉琳， 熊盛青

1 長(zhǎng)安大學(xué)重磁方法技術(shù)研究所， 長(zhǎng)安大學(xué)地質(zhì)工程與測(cè)繪學(xué)院， 長(zhǎng)安大學(xué)西部礦產(chǎn)資源與地質(zhì)工程教育部重點(diǎn)實(shí)驗(yàn)室， 西安 710054 2 中國(guó)自然資源航空物探遙感中心, 北京 100083 3 自然資源部航空地球物理與遙感地質(zhì)重點(diǎn)實(shí)驗(yàn)室, 北京 100083

0 引言

重、磁勘探以其經(jīng)濟(jì)快速、綠色環(huán)保，覆蓋范圍廣、橫向分辨率高等特點(diǎn)在解決地質(zhì)問(wèn)題中得到了廣泛應(yīng)用.為了更好地解決某一地質(zhì)問(wèn)題，往往在不同時(shí)期采集了不同尺度(比例尺、精度)、不同“維度”(衛(wèi)星、航空、地面、海洋及井中等)的重、磁位場(chǎng)數(shù)據(jù).然而，由于存在觀測(cè)精度和觀測(cè)位置及數(shù)據(jù)基準(zhǔn)等差異，同一研究區(qū)采集的位場(chǎng)數(shù)據(jù)存在著精度差異、高度差異及數(shù)據(jù)基準(zhǔn)差異等問(wèn)題.如何將這些不同精度、不同比例尺、不同“維度”和不同基準(zhǔn)的重、磁位場(chǎng)數(shù)據(jù)融合成為統(tǒng)一基準(zhǔn)、同一“維度”(后簡(jiǎn)寫為維度)的位場(chǎng)數(shù)據(jù)就成為重、磁位場(chǎng)數(shù)據(jù)融合的核心問(wèn)題之一.根據(jù)位場(chǎng)數(shù)據(jù)融合中數(shù)據(jù)觀測(cè)比例尺和觀測(cè)維度的差異，本文將位場(chǎng)數(shù)據(jù)融合分為單尺度和多尺度、單維和多維融合問(wèn)題.單尺度指相同網(wǎng)格距、相同比例尺和相同精度的位場(chǎng)數(shù)據(jù)，多尺度指不同網(wǎng)格距、不同比例尺或不同精度的位場(chǎng)數(shù)據(jù)；單維指同一個(gè)連續(xù)觀測(cè)面的位場(chǎng)數(shù)據(jù)，多維指不同高度連續(xù)觀測(cè)面或間斷觀測(cè)面的位場(chǎng)數(shù)據(jù).

單尺度重、磁位場(chǎng)數(shù)據(jù)具有相同的觀測(cè)位置和觀測(cè)精度，不同來(lái)源的單尺度位場(chǎng)數(shù)僅存在數(shù)據(jù)基準(zhǔn)的差異，其融合主要考慮數(shù)據(jù)基準(zhǔn)的統(tǒng)一；而多尺度數(shù)據(jù)由于存在觀測(cè)精度或觀測(cè)比例尺的差異，需要在單尺度數(shù)據(jù)基準(zhǔn)統(tǒng)一的基礎(chǔ)上考慮精度校正.

單尺度數(shù)據(jù)融合方法主要有加權(quán)平均和主成分分析法等加權(quán)類方法(薛重生等，1997；張麗莉等，2003)，其在沒(méi)有基準(zhǔn)差異的單尺度數(shù)據(jù)融合中計(jì)算效果良好，若數(shù)據(jù)基準(zhǔn)存在較大差異，不考慮背景差校正的加權(quán)融合則不滿足位場(chǎng)數(shù)據(jù)融合要求.多尺度數(shù)據(jù)融合方法主要分為兩大類，第一類是插值方法，不同精度的數(shù)據(jù)同時(shí)進(jìn)行插值計(jì)算，主要方法有Kriging插值、加權(quán)最小曲率插值、Shepard插值、樣條插值等(Beylkin and Cramer，2002；劉翔等，2016；劉善偉等，2018a；王虎彪等，2018；柯寶貴等，2018；支澳威和陳華根，2019)，此類方法與單尺度數(shù)據(jù)中的加權(quán)計(jì)算方法思想類似，能夠考慮不同精度數(shù)據(jù)的差異，但直接插值計(jì)算結(jié)果會(huì)受到多尺度數(shù)據(jù)基準(zhǔn)差影響，不同精度數(shù)據(jù)邊界處形成梯級(jí)帶；第二類是由不同分辨率的圖像融合方法延伸而來(lái)，包括基于小波變換(Ranchin and Wald，1993；Deighan and Watts，1997；郭達(dá)志和方濤，1997；Kuroishi and Keller，2005；汪海洪，2005；李雙成等，2006；Douma and De Hoop，2007；Herrmann et al.，2008；吳健生和劉苗，2008；王海江，2012；楊揚(yáng)，2013；畢奔騰等，2016；譚沛森，2018；劉善偉等，2018b；陳勐韜等，2020)和Curvelet變換(Candes and Donoho，2005；Nencinia et al.，2007；張恒磊和劉天佑，2011)的融合方法，這類方法將數(shù)據(jù)分解為不同頻段的信息，通過(guò)對(duì)多來(lái)源數(shù)據(jù)頻段的選擇和重新組合完成數(shù)據(jù)背景差的校正和不同精度數(shù)據(jù)的融合，其在不同分辨率的圖像融合中應(yīng)用效果良好，但由于不同精度的位場(chǎng)數(shù)據(jù)都包含全段頻率的信息，進(jìn)行頻率分解再組合的方法會(huì)損失部分高精度數(shù)據(jù)細(xì)節(jié)信息，不符合位場(chǎng)數(shù)據(jù)優(yōu)選高精度數(shù)據(jù)的原則.

單維數(shù)據(jù)不存在觀測(cè)位置差異，其融合按照上述分類直接考慮單尺度或多尺度融合方法，而多維數(shù)據(jù)則需在單尺度融合和多尺度融合的基礎(chǔ)上再考慮觀測(cè)空間位置的差異.在實(shí)際的數(shù)據(jù)處理工作中，有學(xué)者不考慮觀測(cè)位置的差異，將多維數(shù)據(jù)當(dāng)作單維數(shù)據(jù)進(jìn)行處理，這不符合位場(chǎng)數(shù)據(jù)隨空間高度位置變化的規(guī)律.

多維數(shù)據(jù)融合可以通過(guò)位場(chǎng)延拓方法實(shí)現(xiàn)維度統(tǒng)一后，再進(jìn)行上述單尺度和多尺度數(shù)據(jù)的融合計(jì)算(吳招才等，2018；余春平，2019)；而目前針對(duì)多維數(shù)據(jù)的融合方法主要分為統(tǒng)計(jì)法和解析法兩大類.統(tǒng)計(jì)法主要有最小二乘配置法，利用多維數(shù)據(jù)建立協(xié)方差函數(shù)后進(jìn)行系統(tǒng)偏差估計(jì)和精度估計(jì)，并在此基礎(chǔ)上逐漸發(fā)展了一些迭代算法、正則化方法等(Tscherning and Rapp，1974；Forsberg，1987；楊元喜，1996；Tscherning et al.，1998；Strykowski and Forsberg，1998；Tziavos et al.，1998；Olesen et al.，2002；Kern et al.，2003；Huang et al.，2006；翟振和和孫中苗，2010；劉曉剛等，2012；黃謨濤等，2013b；潘宗鵬等，2013；孫文等，2016；王灼華等，2017；趙啟龍，2017)，該方法能夠?qū)⒉煌^測(cè)面上的多種類型數(shù)據(jù)融合起來(lái)，但需以足夠分辨率的觀測(cè)數(shù)據(jù)為基礎(chǔ)，且經(jīng)驗(yàn)協(xié)方差函數(shù)的選取復(fù)雜、大型矩陣的求解過(guò)程使得該方法適用性和穩(wěn)定性受到限制.許多學(xué)者致力于解析法的研究(Lehmann，1993；Sideris，1996；Li and Sideris，1997；寧津生等，1998；Kern et al.，2003；郝燕玲等，2007；Eicker，2008；Wittwer，2010；黃謨濤等，2013a，2015；王燚等，2015；吳懌昊等，2016；劉金釗等，2020)，該方法在迭代中利用不同觀測(cè)面上解析計(jì)算的位場(chǎng)數(shù)據(jù)殘差修正觀測(cè)數(shù)據(jù)或模型系數(shù)，引入分步平差、推估內(nèi)插、迭代下延、點(diǎn)質(zhì)量等技術(shù)，在大區(qū)域的重、磁力場(chǎng)計(jì)算中應(yīng)用廣泛，能夠省去統(tǒng)計(jì)法中求解大型矩陣逆的問(wèn)題，在迭代解析計(jì)算中完成多維多尺度位場(chǎng)數(shù)據(jù)的融合，但迭代中修正模型數(shù)據(jù)的方法無(wú)法全部保留高精度數(shù)據(jù)，且仍存在迭代穩(wěn)定性的問(wèn)題.

針對(duì)上述研究現(xiàn)狀，本文提出了一套適用于多維多尺度重、磁位場(chǎng)數(shù)據(jù)融合的實(shí)用性方法；對(duì)數(shù)據(jù)維度、基準(zhǔn)、尺度的差異，分別采用位場(chǎng)延拓、回歸分析、精度加權(quán)平均等方法來(lái)解決高度不同、基準(zhǔn)不同、精度不同或維度不同、尺度不同的重、磁位場(chǎng)數(shù)據(jù)融合問(wèn)題，以期得到統(tǒng)一維度、基準(zhǔn)及尺度的重、磁位場(chǎng)數(shù)據(jù)，為進(jìn)一步處理、反演和解釋奠定數(shù)據(jù)基礎(chǔ).

1 多維多尺度重、磁位場(chǎng)數(shù)據(jù)融合方法

為了解決多維多尺度位場(chǎng)數(shù)據(jù)融合問(wèn)題，需研究尺度統(tǒng)一方法、基準(zhǔn)差統(tǒng)一方法(單尺度、多尺度位場(chǎng)數(shù)據(jù)融合方法)和維度統(tǒng)一方法(多維單尺度位場(chǎng)數(shù)據(jù)融合方法).

單尺度數(shù)據(jù)融合方法，采用相關(guān)分析和回歸分析校正基準(zhǔn)差解決數(shù)據(jù)基準(zhǔn)差異；多尺度數(shù)據(jù)融合方法，采用不同精度數(shù)據(jù)加權(quán)解決觀測(cè)尺度和觀測(cè)精度的差異；多維多尺度位場(chǎng)數(shù)據(jù)融合方法，在以上兩類數(shù)據(jù)融合的基礎(chǔ)上采用位場(chǎng)延拓方法解決觀測(cè)維度的差異(圖1).

圖1 多維多尺度重、磁位場(chǎng)數(shù)據(jù)融合計(jì)算流程示意圖Fig.1 The schematic diagram of multi-dimension and multi-scale gravity and magnetic potential data fusion calculation process

1.1 維度統(tǒng)一方法

本文采用位場(chǎng)延拓實(shí)現(xiàn)多維多尺度重、磁位場(chǎng)數(shù)據(jù)融合中的維度統(tǒng)一.本次所使用的位場(chǎng)延拓方法以空間域積分-迭代法(徐世浙，2006)為基礎(chǔ)算法，將不穩(wěn)定的向下延拓問(wèn)題轉(zhuǎn)換為迭代的穩(wěn)定向上延拓問(wèn)題，每次迭代中向上延拓和偏差疊加更新過(guò)程是影響計(jì)算精度和效率的關(guān)鍵，這兩步驟也影響維度統(tǒng)一過(guò)程的計(jì)算精度和效率，因此成為進(jìn)行多維數(shù)據(jù)融合的基礎(chǔ).基于此，本文采用滑動(dòng)窗口加曲面切片技術(shù)(劉芬和王萬(wàn)銀，2019)提高計(jì)算效率，采用向上延拓濾波技術(shù)提高計(jì)算穩(wěn)定性.

(1)位場(chǎng)向上延拓

空間域由原始點(diǎn)(ξ,η,ζ)延拓至計(jì)算點(diǎn)(x,y,z)的基本公式(曾華霖，2005)為：

U(x,y,z)=

(1)

以上解析積分在程序?qū)崿F(xiàn)中使用的是數(shù)值積分的矩形積分公式，則(1)式變?yōu)椋?/p>

(2)

其中，Δξ和Δη分別為x和y方向的點(diǎn)距，H(x,y,z,ξi,ηj,ζij)為核函數(shù)，其值如(3)式，

H(x,y,z,ξi,ηj,ζij)=

本文采用滑動(dòng)窗口和曲面切片的向上延拓計(jì)算技術(shù)(劉芬，2019)提高計(jì)算效率.根據(jù)位場(chǎng)延拓計(jì)算式(式(1))看出，距離越遠(yuǎn)的場(chǎng)源對(duì)計(jì)算點(diǎn)產(chǎn)生的異常影響越小.因此根據(jù)延拓高度在計(jì)算點(diǎn)一定窗口范圍內(nèi)計(jì)算，一般選擇20倍的延拓高度以上，且小于剖面長(zhǎng)度的一半(劉芬等，2019)，可在保證精度的前提下減少計(jì)算量.

對(duì)于延拓高度相同的點(diǎn)(圖2中P點(diǎn)和Q點(diǎn))，在相同的窗口大小內(nèi)，其核函數(shù)矩陣是相等的；平面到平面延拓(圖2由ΓA到ΓS)，可以計(jì)算得到窗口內(nèi)核函數(shù)矩陣，然后直接代入式(2)完成每個(gè)計(jì)算點(diǎn)的延拓.

圖2 向上延拓中滑動(dòng)窗口和曲面切片示意圖Fig.2 The schematic diagram of sliding window and surface slice in upward continuation

若延拓面為曲面(圖2中ΓC)，用一組平面切片Slice將其切割后，采用平面的滑動(dòng)窗口核函數(shù)計(jì)算方法，首先計(jì)算得到每個(gè)切片高度上的固定窗口內(nèi)的核函數(shù)矩陣，落在切片上的計(jì)算點(diǎn)(Q點(diǎn)和M點(diǎn))直接代入，分布在切片之間的計(jì)算點(diǎn)(N點(diǎn))的位場(chǎng)值UN則采用反距離加權(quán)插值得到，其加權(quán)公式(劉芬，2019)為：

(4)

其中，UM、UQ分別為M點(diǎn)和Q點(diǎn)的位場(chǎng)值，hMN、hQN為N點(diǎn)分別與M點(diǎn)、Q點(diǎn)在垂直方向的距離.

由此，通過(guò)在空間域向上延拓計(jì)算中加入滑動(dòng)窗口和曲面切片技術(shù)可在保證精度的前提下大大減少計(jì)算量，提高多維重、磁位場(chǎng)數(shù)據(jù)融合計(jì)算的效率.

設(shè)計(jì)理論模型檢驗(yàn)向上延拓計(jì)算方法技術(shù)的正確性和有效性.場(chǎng)源由五個(gè)直立六面體組成，空間位置和物性參數(shù)(表1)及其水平展布位置(圖3a)已知；模型基礎(chǔ)觀測(cè)曲面高度起伏為-1000～-2500 m(觀測(cè)曲面一，z坐標(biāo)向下為正，圖3b)，觀測(cè)平面為z=0 m，x方向坐標(biāo)范圍0～26000 m，y方向坐標(biāo)范圍0～26000 m，兩個(gè)方向網(wǎng)格間距都為100 m，在觀測(cè)曲面一和z=0 m平面上正演計(jì)算得到重力異常和磁力異常(圖3c、d、e、f).

本次測(cè)試z=0 m平面數(shù)據(jù)向上延拓至觀測(cè)曲面一(-1000～-2500 m)的計(jì)算精度和效率.分別進(jìn)行無(wú)窗口的全區(qū)域延拓計(jì)算、加滑動(dòng)窗口的延拓計(jì)算以及滑動(dòng)窗口結(jié)合曲面切片的延拓計(jì)算，分別統(tǒng)計(jì)了延拓計(jì)算至觀測(cè)曲面一的重力異常(圖4a、b、c)和磁力異常(圖4d、e、f)與理論重、磁力異常(圖3c、d)的均方根誤差、相對(duì)均方根誤差以及計(jì)算時(shí)間(表2，配置Intel(R) Core(i5) CPU 2.80 GHz，內(nèi)存12.0 G的計(jì)算機(jī)).由表2可以看出延拓計(jì)算的重、磁力異常與理論值的偏差都較小，滑動(dòng)窗口和曲面切片的計(jì)算技術(shù)能夠保證延拓計(jì)算的精度；無(wú)窗口的全區(qū)域的延拓計(jì)算耗時(shí)最長(zhǎng)，加滑動(dòng)窗口的延拓計(jì)算耗時(shí)次之，滑動(dòng)窗口結(jié)合曲面切片的延拓計(jì)算耗時(shí)最短，驗(yàn)證了本次采用的滑動(dòng)窗口結(jié)合曲面切片計(jì)算技術(shù)能夠提高計(jì)算效率.這為迭代法向下延拓的多維重、磁位場(chǎng)數(shù)據(jù)融合中維度統(tǒng)一處理奠定了計(jì)算基礎(chǔ).

表1 理論模型體空間位置坐標(biāo)及物性參數(shù)Table 1 The space position coordinates and physical property parameters of theoretical models

表2 理論模型向上延拓計(jì)算誤差統(tǒng)計(jì)表Table 2 The error table of upward continuation calculation of theoretical models

圖3 理論模型平面位置、觀測(cè)面一起伏及重、磁力異常等值線圖(a) 理論模型空間位置平面圖； (b) 觀測(cè)曲面一高程等值線圖； (c) 模型體在觀測(cè)曲面一上引起的重力異常等值線圖； (d) 模型體在觀測(cè)曲面一上引起的磁力異常等值線圖； (e) 模型體在z=0 m平面上引起的重力異常等值線圖； (f) 模型體在z=0 m平面上引起的磁力異常等值線圖.圖中黑色框線表示五個(gè)模型體的邊界位置.Fig.3 The plane position of theoretical models, observation surface One and its gravity and magnetic anomaly contour maps(a) The plane position of theoretical models; (b) The contour map of observation surface One elevation; (c) The contour map of gravity anomaly caused by the models on the observation surface One; (d) The contour map of magnetic anomaly caused by the models on the observation surface One; (e) The contour map of gravity anomaly caused by the models on the z=0 m plane; (f) The contour map of magnetic anomaly caused by the models on the z=0 m plane. The black solid lines represent the edge location of the five models.

圖5 積分-迭代法向下延拓示意圖Fig.5 The schematic diagram of integral-iterative method for downward continuation

(2)位場(chǎng)向下延拓

本文采用積分-迭代法(徐世浙，2006)，將向下延拓轉(zhuǎn)換為迭代法的向上延拓.如圖5所示，將觀測(cè)面ΓB上的位場(chǎng)值UB(x,y,z)向下延拓至ΓA主要進(jìn)行以下幾個(gè)步驟：

(5)

(6)

ΔU(n)′B(x′,y′,z′)=

(7)

其中濾波器H就是向上延拓的核函數(shù)(式(3))；相應(yīng)地，迭代中的校正更新公式(式(6))變?yōu)椋?/p>

圖6 觀測(cè)曲面二地形起伏Fig.6 The contour map of observation surface Two elevation

采用理論模型測(cè)試向下延拓方法技術(shù)的準(zhǔn)確性和精度.觀測(cè)曲面二與觀測(cè)曲面一坐標(biāo)范圍一致，z坐標(biāo)在觀測(cè)曲面一的基礎(chǔ)上抬升3000 m(z坐標(biāo)從-3000 m到-4500 m變化，圖6)，在其正演的重、磁力異常(圖7a、b)基礎(chǔ)上加入各自數(shù)據(jù)幅值5%的隨機(jī)噪聲，即重力異常數(shù)據(jù)噪聲均值為0 mGal、幅值為-0.049～0.055 mGal，磁力異常數(shù)據(jù)噪聲均值為0 nT、幅值為-27.488～27.888 nT，利用加噪后的數(shù)據(jù)(圖7c、d)向下延拓至不同的平面(z=0 m、z=-1000 m、z=-2000 m、z=-3000 m)，得到的向下延拓結(jié)果(圖8)，計(jì)算的均方根誤差和相對(duì)均方根誤差見表3.由圖8可以看出，向下延拓計(jì)算結(jié)果符合位場(chǎng)值隨空間變化的規(guī)律，隨著觀測(cè)面高度與模型體距離的增加，重、磁力異常的幅值減小.由表3可以看出，原始數(shù)據(jù)存在噪聲時(shí)，進(jìn)行低通濾波處理后向下延拓計(jì)算效果良好，各個(gè)延拓結(jié)果的計(jì)算誤差??；且隨著延拓距離的增加，向下延拓計(jì)算精度逐漸降低，這與位場(chǎng)隨空間變化規(guī)律相符合.

圖7 理論模型觀測(cè)曲面二重、磁力異常等值線圖(a) 重力異常； (b) 磁力異常； (c) 加噪聲后的重力異常； (d) 加噪聲后的磁力異常.Fig.7 The contour maps of gravity and magnetic anomaly of theoretical models on the observation surface Two(a) The gravity anomaly; (b) The magnetic anomaly; (c) The gravity anomaly that added noise; (d) The magnetic anomaly that added noise.

圖8 觀測(cè)曲面二向下延拓計(jì)算重、磁力異常等值線圖(a) z=0 m重力異常； (b) z=-1000 m重力異常； (c) z=-2000 m重力異常； (d) z=-3000 m重力異常； (e) z=0 m磁力異常； (f) z=-1000 m磁力異常； (g) z=-2000 m磁力異常； (h) z=-3000 m磁力異常.Fig.8 The contour map of downward continuation results from observation surface Two(a) The gravity anomaly on the z=0 m plane; (b) The gravity anomaly on the z=-1000 m plane; (c) The gravity anomaly on the z=-2000 m plane; (d) The gravity anomaly on the z=-3000 m plane; (e) The magnetic anomaly on the z=0 m plane; (f) The magnetic anomaly on the z=-1000 m plane; (g) The magnetic anomaly on the z=-2000 m plane; (h) The magnetic anomaly on the z=-3000 m plane.

表3 理論模型重、磁力異常向下延拓計(jì)算誤差統(tǒng)計(jì)表Table 3 The error table of gravity or magnetic anomaly downward continuation calculation of theoretical models

同時(shí)，為檢驗(yàn)上述位場(chǎng)延拓計(jì)算方法在多維重、磁位場(chǎng)數(shù)據(jù)融合中的應(yīng)用效果，本文將向下延拓的計(jì)算結(jié)果(圖8)分別與原始理論重、磁力異常進(jìn)行融合(圖9a、c)，同時(shí)對(duì)比不進(jìn)行位場(chǎng)延拓直接“融合”的結(jié)果(即將多維數(shù)據(jù)當(dāng)作單維數(shù)據(jù)直接處理，圖9b、d)，并統(tǒng)計(jì)融合計(jì)算誤差(表4).結(jié)果顯示，在多維數(shù)據(jù)融合計(jì)算中，不考慮數(shù)據(jù)之間的維度差異直接進(jìn)行“融合”是不正確，而本文采用的位場(chǎng)延拓方法可以有效統(tǒng)一不同數(shù)據(jù)的維度差異，計(jì)算誤差小，是進(jìn)行多維多尺度位場(chǎng)數(shù)據(jù)融合的維度統(tǒng)一的有效方法.

表4 理論模型多維重、磁力異常融合計(jì)算誤差統(tǒng)計(jì)表Table 4 The error table of multi-dimension gravity or magnetic anomaly fusion of theoretical models

圖9 多維數(shù)據(jù)融合計(jì)算重、磁力異常立體圖(a) 多維重力異常進(jìn)行維度統(tǒng)一后的融合結(jié)果； (b) 多維重力異常當(dāng)作單維直接進(jìn)行處理的計(jì)算結(jié)果； (c) 多維磁力異常進(jìn)行維度統(tǒng)一后的融合結(jié)果； (d) 多維磁力異常當(dāng)作單維直接進(jìn)行處理的計(jì)算結(jié)果.Fig.9 The 3D surface map of multi-dimension gravity or magnetic anomalies fusion results(a) The fusion results of multi-dimension gravity anomalies after dimension unification; (b) The calculation results of multi-dimension gravity anomalies regarded as single-dimension data; (c) The fusion results of multi-dimension magnetic anomalies after dimension unification; (d) The calculation results of multi-dimension magnetic anomalies regarded as single-dimension data.

1.2 基準(zhǔn)差統(tǒng)一方法

本文采用回歸分析校正基準(zhǔn)的方法技術(shù)解決單尺度數(shù)據(jù)和多尺度數(shù)據(jù)的基準(zhǔn)差異問(wèn)題.不同組數(shù)據(jù)間進(jìn)行回歸分析的前提是其存在顯著性相關(guān)，故需對(duì)重疊區(qū)的多來(lái)源數(shù)據(jù)進(jìn)行相關(guān)顯著性檢驗(yàn).

在單尺度位場(chǎng)數(shù)據(jù)融合中，由于不存在精度或比例尺的差異，在理論上選取任一組數(shù)據(jù)作為基準(zhǔn)數(shù)據(jù)都可，其余數(shù)據(jù)都向此基準(zhǔn)數(shù)據(jù)校正；而對(duì)于存在精度或比例尺差異的多尺度數(shù)據(jù)，須選擇精度最高或比例尺最大的一組數(shù)據(jù)作為校正的基準(zhǔn).

統(tǒng)計(jì)學(xué)領(lǐng)域中，相關(guān)分析是用來(lái)考查兩組地位平等且無(wú)因果關(guān)系的變量相互變化的關(guān)聯(lián)關(guān)系；回歸分析是研究一個(gè)因變量與一個(gè)或多個(gè)自變量之間數(shù)量關(guān)系的統(tǒng)計(jì)方法.應(yīng)用至位場(chǎng)數(shù)據(jù)融合中，相關(guān)分析中為兩組重、磁位場(chǎng)數(shù)據(jù)，例如不同時(shí)期采集的同一地區(qū)的維度相同、尺度相同但存在數(shù)據(jù)基準(zhǔn)差異的重力異常，記為gA、gB，采用Pearson相關(guān)計(jì)算公式，即(9)式，通過(guò)計(jì)算相關(guān)系數(shù)r判定兩組數(shù)據(jù)之間的線性相關(guān)程度(Pearson,1895)：

(9)

不同來(lái)源的重、磁位場(chǎng)數(shù)據(jù)存在的數(shù)據(jù)基準(zhǔn)差主要源自觀測(cè)中基點(diǎn)不同，因此在建立回歸模型過(guò)程中考慮構(gòu)建較簡(jiǎn)單的回歸模型，如一元線性回歸或一元二次回歸即可擬合不同數(shù)據(jù)基準(zhǔn)的差異.求解方法是在建立合適的回歸模型條件下，以最小二乘法求解模型統(tǒng)計(jì)數(shù)，獲得回歸方程，從而依據(jù)其中一組作為基準(zhǔn)的數(shù)據(jù)對(duì)其他位場(chǎng)數(shù)據(jù)進(jìn)行校正.

例如，對(duì)于N元線性回歸，簡(jiǎn)單回歸模型表示為(Galton，1886)：

gA=β0+β1gB+…+βNgB+ε,

(10)

模型中，gA是gB的函數(shù)部分加誤差項(xiàng)，函數(shù)部分反映了數(shù)據(jù)gB與gA之間的變化對(duì)應(yīng)關(guān)系，誤差項(xiàng)ε為隨機(jī)變量，反映了不同時(shí)期采集的位場(chǎng)數(shù)據(jù)gB和gA之間的線性關(guān)系之外的隨機(jī)干擾，β0、β1、…、βN為回歸模型的參數(shù)，由最小二乘估計(jì)得到.

(11)

(12)

根據(jù)這一基本原理，進(jìn)行多組數(shù)據(jù)之間的一元線性回歸分析的處理是選擇一個(gè)回歸的數(shù)據(jù)體gA(x,y)，依次用其他數(shù)據(jù)gB(x,y),gC(x,y),…,gM-1(x,y)與其進(jìn)行回歸校正，并建立簡(jiǎn)單線性回歸模型形式為：

(13)

即只用數(shù)據(jù)本身進(jìn)行校正，得到每個(gè)數(shù)據(jù)體與選取的基準(zhǔn)數(shù)據(jù)的回歸系數(shù)后，進(jìn)行相應(yīng)的計(jì)算即完成了數(shù)據(jù)基準(zhǔn)的回歸分析校正.

圖10 觀測(cè)曲面三地形起伏Fig.10 The contour map of observation surface Three elevation

選取單尺度理論模型數(shù)據(jù)測(cè)試該回歸分析方法在消除基準(zhǔn)差中的正確性和有效性.理論模型場(chǎng)源幾何參數(shù)與物性參數(shù)見表1，為模擬不同時(shí)期采集的存在數(shù)據(jù)基準(zhǔn)差異的不同范圍的單尺度重力、磁力異常數(shù)據(jù)，設(shè)計(jì)兩個(gè)不同觀測(cè)曲面(觀測(cè)曲面一和觀測(cè)曲面三)，觀測(cè)曲面三(圖10)是在觀測(cè)曲面一(圖3b)的基礎(chǔ)上存在空白區(qū)，兩觀測(cè)曲面的重復(fù)區(qū)占總范圍的37%，設(shè)計(jì)觀測(cè)面三原始觀測(cè)異常(數(shù)據(jù)體2)和觀測(cè)面一原始觀測(cè)異常(數(shù)據(jù)體1)之間存在線性的基準(zhǔn)差，其關(guān)系式如式(14)所示：

(14)

對(duì)兩組重力異常(圖11a、b)和磁力異常(圖11c、d)進(jìn)行回歸分析的基準(zhǔn)校正，分別得到校正后重力異常和磁力異常(圖12a、b)和計(jì)算誤差(圖12c、d)，可以看到其校正后重、磁力異常形態(tài)和變化趨勢(shì)符合位場(chǎng)數(shù)據(jù)的特點(diǎn)，且形態(tài)和幅值均與理論重力異常(圖3c)和磁力異常(圖3d)吻合很好，計(jì)算誤差??；重力異常校正計(jì)算的均方根誤差小于2.008×10-4mGal，相對(duì)均方根誤差為0.0079%，磁力異常校正計(jì)算的均方根誤差小于0.135 nT，相對(duì)均方根誤差為0.0087%，計(jì)算精度高，說(shuō)明本次采用的回歸分析能夠?qū)纬叨任粓?chǎng)數(shù)據(jù)基準(zhǔn)差進(jìn)行很好的校正.

多尺度理論模型場(chǎng)源幾何參數(shù)與物性參數(shù)見表1.為模擬多尺度數(shù)據(jù)特點(diǎn)，在觀測(cè)曲面一的基礎(chǔ)上，改變數(shù)據(jù)間隔為50 m并設(shè)置空白區(qū)形成觀測(cè)曲面四(圖13a)，改變數(shù)據(jù)間隔為200 m形成觀測(cè)曲面五(圖13e)，數(shù)據(jù)基準(zhǔn)差關(guān)系式如式(14)，同時(shí)加入不同大小的隨機(jī)噪聲模擬兩者觀測(cè)精度的差異，重力異常數(shù)據(jù)噪聲均值為0 mGal、幅值為-0.028～0.027 mGal(圖13b)和-0.048～0.045 mGal(圖13c)，磁力異常數(shù)據(jù)噪聲均值為0 nT、幅值為-5.826～5.944 nT(圖13f)和-11.861～12.982 nT(圖13g)；兩組數(shù)據(jù)覆蓋范圍不同，存在數(shù)據(jù)重復(fù)區(qū)(數(shù)據(jù)重復(fù)區(qū)占總區(qū)域37%).

對(duì)重力異常、磁力異常分別進(jìn)行回歸分析校正得到原觀測(cè)面上背景差校正后的重、磁力異常(圖13d、h)，可以看到其與理論重力異常(圖3c)和磁力異常(圖3d)的形態(tài)和幅值均吻合很好，基準(zhǔn)校正的效果受觀測(cè)比例尺和隨機(jī)誤差的影響小，重力異常計(jì)算的均方根誤差為0.013 mGal、相對(duì)均方根誤差為0.433%，磁力異常計(jì)算的均方根誤差為3.336 nT、相對(duì)均方根誤差為0.175%，計(jì)算精度較高.同時(shí)，進(jìn)行了不同空白范圍的測(cè)試，即保持觀測(cè)面五的范圍不變，觀測(cè)面四范圍占總區(qū)域比例為9.46%、25.00%、37.86%、47.92%、71.59%，并進(jìn)行誤差統(tǒng)計(jì)(表5)，可以看出兩組數(shù)據(jù)重疊區(qū)不同時(shí)，多尺度位場(chǎng)數(shù)據(jù)基準(zhǔn)校正計(jì)算結(jié)果穩(wěn)定性好，重疊區(qū)覆蓋范圍越廣，校正計(jì)算精度越高.

1.3 尺度統(tǒng)一方法

尺度統(tǒng)一需考慮數(shù)據(jù)觀測(cè)精度或比例尺的不同，在重疊區(qū)需根據(jù)精度或比例尺加權(quán)計(jì)算；為滿足后期數(shù)據(jù)處理和解釋的數(shù)據(jù)需求以及高精度數(shù)據(jù)的完全保留，需將多尺度數(shù)據(jù)統(tǒng)一至最高精度的位場(chǎng)數(shù)據(jù)尺度，本文采用最小曲率方法進(jìn)行網(wǎng)格化(王萬(wàn)銀和邱之云，2011).

重、磁位場(chǎng)數(shù)據(jù)是地下密度、磁性不均勻體的綜合反映，高精度數(shù)據(jù)包含區(qū)域場(chǎng)信息和更豐富的細(xì)節(jié)信息，為在數(shù)據(jù)融合能保留高精度數(shù)據(jù)，在數(shù)據(jù)重疊區(qū)精度加權(quán)的原則為：高精度數(shù)據(jù)權(quán)重為1，低精度權(quán)重為0；若精度相同，則權(quán)重相同(即取平均值).需要注意的是，在進(jìn)行位場(chǎng)延拓時(shí)，延拓的距離越大，計(jì)算的誤差相應(yīng)也越大，因此不同高度數(shù)據(jù)經(jīng)位場(chǎng)延拓?fù)Q算至同一個(gè)平面或曲面后可信度也有所差別，所以，多維單尺度數(shù)據(jù)重疊區(qū)也存在精度差異，需要加權(quán)處理.

圖11 單尺度原始重力、磁力數(shù)據(jù)等值線圖(a) 重力異常一； (b) 重力異常二； (c) 磁力異常一； (d) 磁力異常二.白色區(qū)域?yàn)閿?shù)據(jù)空白區(qū).Fig.11 The contour map of original single-scale gravity and magnetic data(a) The first gravity anomaly; (b) The second gravity anomaly; (c) The first magnetic anomaly; (d) The second magnetic anomaly. The white areas indicate that there has no data in the areas.

圖12 單尺度重力、磁力異?；貧w分析計(jì)算結(jié)果及誤差等值線圖(a) 重力異?；貧w分析結(jié)果； (b) 磁力異常回歸分析結(jié)果； (c) 重力異常回歸分析計(jì)算誤差； (d) 磁力異常回歸分析計(jì)算誤差.Fig.12 The contour map of regression results of single-scale gravity and magnetic anomalies(a) The regression results of gravity anomalies; (b) The regression results of magnetic anomalies; (c) The error of regression results of gravity anomalies; (d) The error of regression results of magnetic anomalies.

圖13 多尺度觀測(cè)面、原始數(shù)據(jù)及回歸分析結(jié)果等值線圖(a) 觀測(cè)面四； (b) 觀測(cè)面四重力異常； (c) 觀測(cè)面五重力異常； (d) 重力異常校正結(jié)果； (e) 觀測(cè)面五； (f) 觀測(cè)面四磁力異常； (g) 觀測(cè)面五磁力異常； (h) 磁力異常校正結(jié)果.白色區(qū)域?yàn)閿?shù)據(jù)空白區(qū).Fig.13 The contour map of multi-scale observation surface, raw data and regression result(a) The observation surface Four; (b) The gravity anomaly on observation surface Four; (c) The gravity anomaly on observation surface Five; (d) The regression results of gravity anomalies; (e) The observation surface Five; (f) The magnetic anomaly on observation surface Four; (g) The magnetic anomaly on observation surface Five; (h) The regression results of magnetic anomalies. The white areas indicate that there has no data in the areas.

表5 理論模型多尺度重、磁力異?；貧w分析校正誤差統(tǒng)計(jì)表Table 5 The error table of multi-scale gravity or magnetic anomaly regression of theoretical models

根據(jù)多尺度位場(chǎng)數(shù)據(jù)融合的方法原理，選取已經(jīng)完成維度統(tǒng)一和基準(zhǔn)差統(tǒng)一的單尺度數(shù)據(jù)(圖12)和多尺度數(shù)據(jù)(圖13d、h)，進(jìn)行單尺度和多尺度融合(圖14).由圖14可以看出，位場(chǎng)數(shù)據(jù)融合計(jì)算結(jié)果能夠保留高精度、大比例尺數(shù)據(jù)的細(xì)節(jié)信息，同時(shí)對(duì)研究區(qū)內(nèi)空白區(qū)進(jìn)行了有效補(bǔ)充，不同數(shù)據(jù)接邊位置光滑，沒(méi)有梯級(jí)帶和突跳.同時(shí)，為測(cè)試重疊區(qū)范圍對(duì)多尺度數(shù)據(jù)融合結(jié)果的影響，設(shè)計(jì)不同重疊區(qū)模型測(cè)試并統(tǒng)計(jì)了融合計(jì)算結(jié)果的誤差(表6)，重疊區(qū)覆蓋范圍不同時(shí)，融合計(jì)算誤差均很小，相對(duì)均方根誤差均小于1%，計(jì)算穩(wěn)定；且重疊區(qū)范圍越廣，融合精度越高.因此，本文提出的多尺度位場(chǎng)數(shù)據(jù)融合方法在數(shù)據(jù)重疊區(qū)范圍不同時(shí)均可得到較好的融合計(jì)算結(jié)果.

最終設(shè)計(jì)多維多尺度理論模型測(cè)試其位場(chǎng)數(shù)據(jù)融合計(jì)算的正確性和有效性.場(chǎng)源參數(shù)如表1，設(shè)計(jì)3個(gè)不同高度、不同尺度、數(shù)據(jù)基準(zhǔn)有差異的不同覆蓋范圍的觀測(cè)面，觀測(cè)面一(圖3b)，觀測(cè)面六(在觀測(cè)面一的基礎(chǔ)上z坐標(biāo)增大1000 m，z坐標(biāo)范圍為0～-1500 m，設(shè)置空白區(qū))，觀測(cè)面七(在觀測(cè)面一的基礎(chǔ)上z坐標(biāo)減小2000 m，z坐標(biāo)范圍為-3000～-4500 m；網(wǎng)格間距調(diào)整為200 m并設(shè)置空白區(qū))，以其理論重力異常(圖15b)和磁力異常(圖15c)進(jìn)行融合測(cè)試.

采用本文的多維多尺度重、磁位場(chǎng)數(shù)據(jù)融合方法.首先將多維數(shù)據(jù)以維度統(tǒng)一方法計(jì)算至同一個(gè)觀測(cè)面上變?yōu)閱尉S數(shù)據(jù)，接著分別統(tǒng)一單維多尺度數(shù)據(jù)的基準(zhǔn)差異和尺度差異，得到該觀測(cè)面上的融合結(jié)果，利用位場(chǎng)延拓得到了z=0 m、z=-2000 m、z=-4000 m和z=-6000 m四個(gè)平面的重、磁力異常(圖16).可以看到其融合結(jié)果與各個(gè)計(jì)算面的理論重、磁力異常均吻合較好，重、磁力異常等值線光滑、異常變化趨勢(shì)一致；隨著延拓高度的增加，異常形態(tài)變化趨于平緩，符合位場(chǎng)數(shù)據(jù)隨高度變化的規(guī)律；融合計(jì)算結(jié)果與理論異常的誤差統(tǒng)計(jì)表(表7)顯示，其計(jì)算至不同高度平面結(jié)果均較穩(wěn)定，計(jì)算精度較高，驗(yàn)證了該方法在多維多尺度重磁位場(chǎng)數(shù)據(jù)融合中的有效性.

圖14 單尺度和多尺度重、磁位場(chǎng)數(shù)據(jù)融合結(jié)果等值線圖(a) 單尺度重力異常融合結(jié)果； (b) 單尺度磁力異常融合結(jié)果； (c) 多尺度重力異常融合結(jié)果； (d) 多尺度磁力異常融合結(jié)果.Fig.14 The contour map of fusion result of single-scale and multi-scale gravity and magnetic anomalies(a) The fusion result of single-scale gravity anomalies; (b) The fusion result of single-scale magnetic anomalies; (c) The fusion result of multi-scale gravity anomalies; (d) The fusion result of multi-scale magnetic anomalies.

表6 理論模型多尺度重、磁力異常融合計(jì)算誤差統(tǒng)計(jì)表Table 6 The error table of multi-scale gravity or magnetic anomaly fusion of theoretical models

圖15 多維多尺度原始觀測(cè)面起伏Fig.15 The 3D surface map of original multi-dimension multi-scale observation surface

圖16 多維多尺度重、磁異常及融合計(jì)算結(jié)果立體圖(a) 觀測(cè)面重力異常； (b) 觀測(cè)面磁力異常； (c) 重力異常融合結(jié)果； (d) 磁力異常融合結(jié)果.Fig.16 The 3D surface map of multi-dimension multi-scale anomalies and fusion results of gravity and magnetic anomalies(a) The gravity anomalies on the observation surface; (b) The magnetic anomalies on the observation surface; (c) The fusion results of gravity anomalies; (d) The fusion results of magnetic anomalies.

表7 理論模型多維多尺度重、磁力異常融合誤差統(tǒng)計(jì)表Table 7 The error table of multi-dimension multi-scale gravity or magnetic anomaly fusion of theoretical models

2 實(shí)際資料處理試驗(yàn)

為檢驗(yàn)本文多維多尺度重、磁位場(chǎng)數(shù)據(jù)融合方法實(shí)用性，選取在中國(guó)某地航空和地面觀測(cè)的不同比例尺重力異常數(shù)據(jù)進(jìn)行多維多尺度位場(chǎng)數(shù)據(jù)融合試驗(yàn).

2.1 單維多尺度重力數(shù)據(jù)融合試驗(yàn)

選取中國(guó)某地1∶50000地面實(shí)測(cè)重力異常(數(shù)據(jù)體一，圖17b)和1∶200000地面實(shí)測(cè)重力異常數(shù)據(jù)(數(shù)據(jù)體二，圖17c)檢驗(yàn)單維多尺度重、磁位場(chǎng)數(shù)據(jù)融合方法的實(shí)用效果.兩個(gè)不同比例尺的地面數(shù)據(jù)的觀測(cè)面均為地形面(圖17a)，地形高度起伏變化范圍是747～2623 m，兩次重力觀測(cè)得到了不同比例尺的數(shù)據(jù)體一(幅值-8～100 mGal)和數(shù)據(jù)體二(幅值-43～96 mGal)，網(wǎng)格化后X和Y方向的網(wǎng)格間距都為150 m，其異常趨勢(shì)和幅值變化在相同位置上大體一致，但1∶50000的數(shù)據(jù)反映更多細(xì)節(jié)信息，而1∶200000數(shù)據(jù)的變化趨勢(shì)相對(duì)光滑.

利用單維多尺度重、磁位場(chǎng)數(shù)據(jù)融合方法將兩組不同比例尺的重力數(shù)據(jù)進(jìn)行融合，以高精度大比例尺數(shù)據(jù)為基準(zhǔn)，將小比例尺數(shù)據(jù)進(jìn)行數(shù)據(jù)基準(zhǔn)差校正，相關(guān)系數(shù)為0.883，回歸模型的斜率為0.996，截距(即數(shù)據(jù)基準(zhǔn)差)為5.636 mGal；校正后的小比例尺數(shù)據(jù)與大比例尺數(shù)據(jù)進(jìn)行偏差統(tǒng)計(jì)，其偏差平均值為0.017 mGal，均方根偏差為7.546 mGal，相對(duì)均方根偏差6.968%，其偏差較小能夠保證校正后數(shù)據(jù)與高精度數(shù)據(jù)的一致性；之后進(jìn)行數(shù)據(jù)融合得到單維多尺度重力數(shù)據(jù)融合結(jié)果(圖17d)，融合結(jié)果重力異常的形態(tài)和幅值(33～100 mGal)均與原觀測(cè)面基本一致，保留了地面高精度的大比例尺數(shù)據(jù)(圖17b)所包含的細(xì)節(jié)信息，同時(shí)用較低精度的小比例尺數(shù)據(jù)有效補(bǔ)充了高精度大比例尺數(shù)據(jù)的空白區(qū)，融合后的重力異常也符合位場(chǎng)的特點(diǎn)，驗(yàn)證了該融合方法在不同尺度實(shí)測(cè)重力異常的融合中的實(shí)用性.

圖17 地形起伏、觀測(cè)重力異常及融合重力異常等值線圖(a) 地形起伏； (b) 數(shù)據(jù)體一1∶50000地面實(shí)測(cè)重力異常； (c) 數(shù)據(jù)體二1∶200000地面實(shí)測(cè)重力異常； (d) 單維多尺度融合重力異常；黑色點(diǎn)代表觀測(cè)點(diǎn).Fig.17 The contour map of topography, measured gravity anomaly and the fusion result(a) The topography; (b) The gravity anomaly measured on the ground at scale 1∶50000; (c) The gravity anomaly measured on the ground at scale 1∶200000; (d) The fusion result of single-dimension multi-scale gravity anomaly; The black points represent the location of the observation points.

2.2 多維多尺度重力數(shù)據(jù)融合試驗(yàn)

為檢驗(yàn)多維多尺度重、磁位場(chǎng)數(shù)據(jù)融合方法的實(shí)用性，在地面不同比例尺重力異常的基礎(chǔ)上增加1∶50000航空實(shí)測(cè)重力異常數(shù)據(jù)，航空觀測(cè)面GPS高度為1638～2705 m(圖18a)，測(cè)線間隔為500 m，測(cè)點(diǎn)間隔為20 m，重力異常幅值為-43～96 mGal(圖18b)，X和Y方向的網(wǎng)格化尺寸為150 m，航空重力異常與地面重力異常變化趨勢(shì)大體一致，但總體幅值較低；兩者比較看出，地面數(shù)據(jù)顯示更多細(xì)節(jié)信息，航空數(shù)據(jù)觀測(cè)比例尺雖與地面1∶50000重力數(shù)據(jù)一致，但其觀測(cè)精度低于地面1∶50000數(shù)據(jù)，航空數(shù)據(jù)主要在覆蓋范圍和工作效率上更具有優(yōu)勢(shì).

采用本文多維多尺度重、磁位場(chǎng)數(shù)據(jù)融合方法對(duì)數(shù)據(jù)體一(1∶50000地面重力異常，圖17b)、數(shù)據(jù)體二(1∶200000地面重力異常，圖17c)和數(shù)據(jù)體三(1∶50000航空重力異常，圖18b)三組航空和地面實(shí)測(cè)的不同精度重力異常數(shù)據(jù)進(jìn)行融合試驗(yàn).

數(shù)據(jù)體一和數(shù)據(jù)體三融合過(guò)程中，航空數(shù)據(jù)采用位場(chǎng)延拓計(jì)算至地面(維度統(tǒng)一)后進(jìn)行基準(zhǔn)差統(tǒng)一，與地面數(shù)據(jù)對(duì)比，其平均偏差小于10-5mGal，均方根偏差為7.379 mGal，相對(duì)均方根偏差為6.813%.可以看出，數(shù)據(jù)體三計(jì)算至地面后的精度略高于數(shù)據(jù)體二(相對(duì)均方根偏差為6.968%)，但精度相差不大；數(shù)據(jù)融合重力異常如圖19a所示，融合重力異常符合位場(chǎng)數(shù)據(jù)特征，形態(tài)和幅值均與原觀測(cè)面基本一致，幅值范圍為-44～124 mGal，保留了地面高精度數(shù)據(jù)(圖17b)包含的細(xì)節(jié)信息，同時(shí)用效率高、覆蓋范圍廣的航空重力數(shù)據(jù)有效地補(bǔ)充了地面重力數(shù)據(jù)的空白區(qū)，兩者數(shù)據(jù)相接的位置趨勢(shì)一致且光滑無(wú)梯級(jí)帶.

根據(jù)上述兩組數(shù)據(jù)融合的精度比較結(jié)果，該地區(qū)三組多維多尺度數(shù)據(jù)融合時(shí)，首先將計(jì)算至地面后精度相差不大的航空1∶50000和地面1∶200000重力異常數(shù)據(jù)融合，加權(quán)平均的權(quán)重為0.5，得到初步融合結(jié)果(圖19b)，再將其與精度最高的地面1∶50000重力數(shù)據(jù)進(jìn)行多尺度融合，即為此三組多維多尺度數(shù)據(jù)融合重力異常(圖19c).

數(shù)據(jù)體二和三的融合結(jié)果(圖19c)已經(jīng)結(jié)合了大范圍的航空重力數(shù)據(jù)與小比例尺小范圍地面數(shù)據(jù)的各自優(yōu)勢(shì)，經(jīng)校正后異常形態(tài)趨勢(shì)一致，最終融合結(jié)果(圖19c)相對(duì)于初步融合結(jié)果又包含更多地面高精度細(xì)節(jié)信息；統(tǒng)計(jì)初步融合結(jié)果校正至地面高精度數(shù)據(jù)后的偏差，平均值小于10-5mGal，均方根偏差為6.556 mGal，相對(duì)均方根偏差為6.054%，由此看出三組數(shù)據(jù)整體融合的效果優(yōu)于任意兩組數(shù)據(jù)，同樣驗(yàn)證了該思路和方法在多維多尺度位場(chǎng)數(shù)據(jù)融合中的有效性.

圖18 航空觀測(cè)GPS高度和實(shí)測(cè)重力異常等值線圖(a) 觀測(cè)面GPS高度； (b) 1∶50000實(shí)測(cè)航空重力異常.Fig.18 The contour map of observed surface GPS altitude and gravity anomaly measured by aviation(a) The GPS altitude of observation surface; (b) Measured gravity anomaly at scale 1∶50000 by aviation.

圖19 多維多尺度融合重力異常等值線圖(a) 航空1∶50000與地面1∶50000重力異常融合； (b) 航空1∶50000與地面1∶200000重力異常融合； (c) 航空1∶50000與地面1∶200000、地面1∶50000重力異常融合； (d) 航空1∶50000與地面1∶200000、地面1∶50000重力異常當(dāng)作單維直接處理結(jié)果.Fig.19 The contour map of fusion results of multi-dimension multi-scale gravity anomaly(a) The fusion result of aviation 1∶50000 and ground 1∶50000 gravity anomaly; (b) The fusion result of aviation 1∶50000 and ground 1∶200000 gravity anomaly; (c) The fusion result of aviation 1∶50000, ground 1∶200000 and ground 1∶50000 gravity anomaly; (d) The processing result of aviation 1∶50000, ground 1∶200000 and ground 1∶50000 gravity anomaly regarded as single-dimension data.

另外，為了進(jìn)一步驗(yàn)證和對(duì)比多維數(shù)據(jù)融合中維度統(tǒng)一的必要性和有效性，本文測(cè)試將多維數(shù)據(jù)當(dāng)作單維數(shù)據(jù)不進(jìn)行位場(chǎng)延拓直接“融合”(圖19d)，其結(jié)果顯示在不同數(shù)據(jù)的接邊處存在明顯梯級(jí)帶，且異常值的變化趨勢(shì)無(wú)法達(dá)到有效一致，不符合位場(chǎng)變化規(guī)律，由此也驗(yàn)證不考慮觀測(cè)位置高度變化而將多維數(shù)據(jù)當(dāng)作單維數(shù)據(jù)直接“融合”的做法是錯(cuò)誤的.本文多維多尺度重、磁位場(chǎng)數(shù)據(jù)融合方法原理及計(jì)算技術(shù)對(duì)于不同觀測(cè)比例尺的地面和航空數(shù)據(jù)進(jìn)行了有效融合，驗(yàn)證了該方法在實(shí)測(cè)地面和航空數(shù)據(jù)融合中的實(shí)用性.

3 結(jié)論與建議

本文研究了多維多尺度重、磁位場(chǎng)數(shù)據(jù)融合問(wèn)題，采用改進(jìn)的積分-迭代位場(chǎng)延拓方法實(shí)現(xiàn)維度統(tǒng)一，采用相關(guān)分析和回歸分析實(shí)現(xiàn)數(shù)據(jù)基準(zhǔn)統(tǒng)一，采用數(shù)據(jù)精度加權(quán)平均的方法實(shí)現(xiàn)尺度統(tǒng)一，解決了多維多尺度重、磁位場(chǎng)數(shù)據(jù)融合中維度、尺度和數(shù)據(jù)基準(zhǔn)差異的問(wèn)題，提出了一套適用于多維多尺度重、磁位場(chǎng)數(shù)據(jù)融合的方法.

理論模型測(cè)試和實(shí)測(cè)資料試驗(yàn)顯示，該方法適用于航空、地面或海洋觀測(cè)的重、磁位場(chǎng)數(shù)據(jù)融合，融合重、磁力異常能夠作為后續(xù)位場(chǎng)數(shù)據(jù)處理、反演和解釋的基礎(chǔ)數(shù)據(jù)，具有良好的實(shí)用價(jià)值.本文提出的重、磁位場(chǎng)數(shù)據(jù)融合的思路也可推廣至重、磁力異常各分量或張量數(shù)據(jù)的融合處理中.

致謝感謝論文評(píng)審專家提出的修改意見和編輯對(duì)論文的編輯、加工.