海洋磁測中單臺站地磁數據補缺方法探討

張向宇 關永賢 王功祥

（1.國土資源部海底礦產資源重點實驗室 廣州海洋地質調查局，廣東 廣州 510075；2.南方海洋科學與工程廣東省實驗室（廣州），廣東 廣州 511458）

隨著國家海洋戰略的實施，海洋測量工作量與日俱增，海洋磁力測量已成為我國涉海部門進行海洋探測的重要手段之一。日變改正是海洋磁測數據處理中一個重要環節，而地磁臺站日變數據的質量直接影響著日變改正的效果，進而影響最終提交的磁測成果數據的精度。一般是選擇工區附近的固定地磁臺站或在大范圍地磁普測期間在工區內投放海底地磁觀測站來獲取地磁日變數據，但受制于測量環境等因素的影響，經常會出現部分地磁日變數據缺失的情況，導致部分測線無法進行日變改正的情況，根據海洋調查規范[1]要求，對于無法進行日變改正的測線需作廢處理，這樣將導致工作量不達標，進而影響成果的提交。

目前國內外學者對于日變數據的補缺處理大多是運用多臺站數據進行計算，單汝儉等[2]在1990年提出了二維多項式最小二乘擬合法、時空擬合法和線性內插法這三種局部地區地磁日變的擬合方法；邊剛等[3]在2009年對加權平均法和函數擬合法進行了分析，并進行了實例驗證，該方法不考慮經度效應的影響，使用緯度方向的距離來計算某點的日變改正值；卞光浪等[4]在2010年利用緯差加權法進行多站地磁日變改正值計算，將日變數據的差異只建立在緯度差關系上，通過實例進行計算驗證，證實該方法較距離加權法效果更好；張向宇等[5]在2016年在緯差加權法的基礎上引入回歸分析法推算某處地磁日變數據，通過實例分析，驗證了方法的可靠性；彭飛等[6]在2015年對調和分析法進行了分析，建立了日變數據處理諧波分析模型，將磁靜和磁擾數據進行了合理分離；劉帆等[7]在2016年對諧波分析法進行了進一步研究，采用最小二乘法確定諧波系數，再取一定截斷階數的傅里葉級數作為太陽靜日變化模型，以此進行磁擾數據的分離，取得了較好效果。

在實際工作中，經常會遇到獲取的某一地磁臺站缺失部分數據的情況，這時會考慮用上述文獻中提到的多種地磁日變數據計算方法進行補缺，但這些方法均要求有合適的地磁臺站數據作為樣本，若待補缺地磁臺站所處位置附近無可用地磁臺站數據作為樣本，則依然無法完成補缺工作，最終會直接影響日變改正工作的進行。因此針對單一地磁臺站缺失部分數據的情況，研究相關方法進行數據補缺，尤其是在附近無可用地磁臺站時進行數據補缺是十分必要的。

1 方法設計

本方法結合了信號分析中的諧波分析法和統計學中的聚類分析和相關分析法，首先使用諧波分析法對磁測數據進行頻譜分離，再使用相關分析和聚類分析對靜日變數據進行遴選和計算，獲取靜日變樣本數據，最后與分離出的磁擾數據合成，進行數據重構，得到最終的結果，從而完成數據補缺，具體的流程如圖1所示。

圖1 單臺站地磁數據補缺方法流程圖

1.1 數據準備

需要準備的數據有兩個，一個是待補缺臺站的數據，即除了缺失數據時間段外的所有數據，此部分數據為必需數據，另一個是參考地磁臺站數據，只需提供參考地磁臺站在缺失數據時間段內的數據，若缺失時間段為磁靜日，則該數據可不必提供，磁靜日的判斷需參考地磁指數，一般kp指數大于3-即為磁擾發生時。參考地磁臺站需與待補缺地磁臺站處于相同緯度，且經度差距越小越好，若相同緯度處無地磁臺站數據可用，則以緯度最近原則，選擇與待補缺臺站距離最近的地磁臺站為參考臺站，參考臺站需在缺失數據時間段內數據完備齊全。

1.2 數據分離

需要對準備好的待補缺臺站數據和參考臺站數據分別進行靜日變化和磁擾的分離，采用的方法是諧波分析，首先對準備好的待補缺臺站已有數據進行分析，選擇已有數據中的磁靜日數據進行數據分離，逐日分離出磁靜數據備用，對參考地磁臺站數據同樣進行數據分離，并保留待補缺時間段內的磁擾數據。

1.3 樣本數據選取

靜日變化數據是以天為單位具有周期性的，可通過諧波分析法對周期性磁靜日地磁變化規律進行重構，但重構時受測量噪音或低幅度磁擾的影響，會導致分離出的靜日變化數據不準確，為了減小這種噪音影響，采用以多天的靜日數據為樣本，求取均值的方法獲取最終的靜日數據，可有效減小誤差，這里在做樣本選擇時運用的方法是相關分析及聚類分析。

1.4 數據合成

對于選擇好的靜日變化樣本數據，以天為單位對其求均值可獲取最終的靜日變化，因靜日變化數據是以天為周期重復的，故只需將獲得的磁靜數據依完整的一天時間為單位進行復制即可，若待補缺時間為磁平靜期，則補缺工作到此完成，若待補缺時間為磁擾期，則需將由參考地磁臺站數據中分離出的磁擾數據依全球時附加到靜日變化數據中，合成數據后，即可完成數據的補缺。

2 實例分析

某工區野外測量時間為2019年7月4日至27日，作業期間投放了兩臺海底地磁觀測站，兩臺海底觀測站投放位置緯度相同，經度相差17度，設定西側的海底地磁觀測站為1號臺站，另一臺為2號臺站，兩臺站實測為F分量，數據間隔均為1 min。假定1號臺站7月21日缺失地磁日變數據，則以2號臺站做為參考臺站，采用上述方法進行數據補缺，并將計算結果與1號臺站21日的實測數據進行對比。

首先需要對1號臺站數據采用諧波分析法進行數據分離，分離時取諧波次數為5，時間間隔為1 min，測量時間內kp和ap指數如圖2所示。

由圖2可知2019年7月至9日及7月14日至20日這9天均為磁平靜日，分離后的磁靜數據如圖3所示。將分離出的9天的磁靜數據進行相關分析，結果如表1所示。

根據相關分析的結果，選擇相關性較高的9日、15日、16日、20日數據為樣本，其中20日最接近21日，因此以20日的靜日變數據為基準，通過聚類分析找尋與20日靜日變數據最接近的數據，對上述日期分離后的靜日變數據進行聚類分析，分析后的樹狀圖如圖4所示[8-10]。

圖2 測量時間內地磁指數分布圖

圖3 1號臺站磁靜數據分離結果對比圖

表1 磁靜數據相關性分析表

圖4 1號臺站靜日變數據聚類分析樹狀圖

從圖4中看出，與20日的靜日變數據聚類最接近的是14日的靜日變數據，且從表1中可以看出14日和20日的磁靜數據相關性也較高，相關指數達到0.946，因此最終選擇14日和20日靜日變數據為最終樣本，對其求取均值后，得到計算出的靜日變數據，如圖5所示。

圖5 1號臺站靜日變化計算結果曲線

通過查詢kp指數，7月21日的kp指數大于3-，為磁擾日，因此需要從2號地磁臺21日的實際觀測數據中分離出磁擾數據，在分離磁擾數據前，先選取7月24日和25日這兩個磁擾日1號臺站及2號臺站的觀測數據進行比較，結果見圖6。

圖6 兩臺站磁擾日地磁日變觀測數據對比圖

從圖6中看出，兩臺站磁擾發生時間一致，無時延，磁擾發生幅度相近，量級一致，對這兩組數據進行相關性分析，求得相關系數為0.938，相關性較好。由此說明2號臺站作為1號臺站磁擾日數據補缺的參考地磁臺站是合理的。對2號臺站7月21日數據進行數據分離，得到磁擾數據，結果見圖7所示。

將由1號臺站數據計算的靜日變數據和由2號臺站分離出的磁擾數據進行組合后，得到1號臺站7月21日的地磁日變數據，即完成補缺工作，結果見圖8。圖中紅色曲線為實際觀測值，藍色曲線為補缺數據。

圖7 2號臺站分離的磁擾曲線

圖8 1號臺站7月21日數據對比圖

統計計算數據和實測數據的偏差，最大值為19.2 nT，均值為8.2 nT，因7月21日是磁擾發生時段，而磁擾的發生本身對計算會帶來很多影響，海洋調查規范[1]中要求磁異常ΔT精度≤5 nT，即原始磁測數據經過預處理及精細處理后得到的測網磁異常ΔT交點差中誤差≤5 nT，日變改正為預處理中的一項，以上計算數據與實測數據偏差均值為8.2 nT，即日變改正誤差為8.2 nT，這種程度的誤差可在后續數據調差處理中進行一定程度的消除，使得最終的磁異常ΔT，值精度達標。另外從圖8中看出，兩條曲線形態具有一致性，但存在時延和幅值差異，這是因為靜日數據是由該臺站其他日期靜日數據合成而來，靜日數據雖是以天為單位的周期函數，但在每個周期內受地磁場變化影響，幅值相位會有所不同，尤其是每日0時至12時段內，是地磁日變化幅度較大的時段，這種偏差會更為明顯，而這種因地磁場變化帶來的幅值和相位的偏差也是不可避免的，在靜日數據合成前采用相關分析也是為了盡量減小這種偏差。通過以上實例分析可知，采用本文方法重構的地磁日變數據與實測日變數據曲線形態一致，兩組數據偏差在合理范圍內，說明該方法在實測數據處理中具有一定的效果，可有效解決單臺站數據缺失的問題。

3 結論

本文通過對方法及實例的分析研究，得到以下結論：（1）該方法結合了信號分析中的諧波分析法和統計學中的聚類分析法，這兩種方法都是較為成熟的方法，易于實現，且可有效提升結果的準確性；（2）當待補缺時間為地磁平靜時，則不需其它臺站數據做參考，便可完成數據補缺；（3）該方法在使用時需要參考地磁指數來區分磁平靜時和磁擾時，再依據地磁日變數據的特性，對數據進行組合，得到最終的結果。

綜上，本文使用的數據補缺方法在可獲取數據有限的情況下，尤其是無距離相近臺站數據可用時，可有效地解決單一地磁臺站數據缺失的問題。在磁平靜時可只憑借臺站本身已有數據便可完成數據補缺，解決了日變改正中時常遇到的問題，保證了磁測數據預處理工作的進行。