眭曉虹 張潤寧 萬曉云 李洋

(1 航天東方紅衛(wèi)星有限公司，北京 100094) (2 錢學森空間技術實驗室，北京 100094) (3 北京空間飛行器總體設計部，北京 100094)

?

基于衛(wèi)星測高數(shù)據(jù)的海底地形反演及分析

眭曉虹1,2張潤寧3萬曉云2李洋1

(1 航天東方紅衛(wèi)星有限公司，北京 100094) (2 錢學森空間技術實驗室，北京 100094) (3 北京空間飛行器總體設計部，北京 100094)

由于受限于觀測條件，海底地形模型在分辨率與精度兩方面較陸地和大地水準面模型有很大的差距。傳統(tǒng)的海底地形測量是利用聲吶船測技術獲取水深數(shù)據(jù)，但該測量方法測量效率低、測量點分布稀疏不均勻。文章詳細研究并推導了基于衛(wèi)星測高數(shù)據(jù)進行海底地形反演的方法，采用衛(wèi)星測高重力異常數(shù)據(jù)和少量船測水深數(shù)據(jù)，對北太平洋部分區(qū)域進行了海底地形反演，得到1′×1′分辨率的海底地形模型，并對反演關系中的非線性項進行了分析。仿真結果驗證了此方法的可行性，指出海底地形與重力異常數(shù)據(jù)在一定波段內(nèi)以線性關系為主，疊加非線性二次項和三次項可較好地提升反演精度，高于三次項的部分對反演的結果沒有明顯的作用，可以忽略。

海底地形；衛(wèi)星測高；重力異常；聲吶船測技術

1 引言

海底地形是描述海洋水下固體表面起伏的物理概念，對海底地形的精確探測不僅對研究地球物理科學具有重大的意義，同時對潛艇水下導航等軍事活動、海底管線鋪設等經(jīng)濟活動都具有指導作用[1-2]。在傳統(tǒng)的海洋測繪中，聲吶船測技術是獲取測量船到海底地形起伏面的深度(即水深)數(shù)據(jù)的核心手段。利用測量船搭載聲吶技術測量水深具有精度高的優(yōu)點，但缺點是效率低、測量點分布稀疏不均勻，無法有效地實現(xiàn)高分辨率、高精度、全球化海底地形的探測。

自20世紀70年代發(fā)射第一顆搭載雷達高度計的測高衛(wèi)星起，應用衛(wèi)星測高技術獲取全球海洋大地水準面高的數(shù)據(jù)變?yōu)榱丝赡躘3-5]。大地水準面作為最接近全球平均海面的重力場等勢面，其形狀反映了地球內(nèi)部物質構造、密度及分布等信息，與海底地形在一定程度上具有高相關性[6]。Dixon、Ramillien、Hwang、黃謨濤等多名國內(nèi)外研究人員基于衛(wèi)星測高數(shù)據(jù)提取海面重力測量數(shù)據(jù)，并對部分海域進了海底地形反演[7-10]，但這些反演過程只考慮海面重力數(shù)據(jù)與海底地形線性相關性，對于非線性的相關程度未予考慮。

本文基于Smith&Sandwell(S&S)的經(jīng)典方法[1]，系統(tǒng)探討了應用衛(wèi)星測高數(shù)據(jù)進行海底地形反演的方法，并進行了仿真分析。仿真結果表明，基于S&S方法可以有效地反演海底地形。針對前人忽略的反演關系中非線性項部分，文中亦做了仿真分析，結果表明，海底地形與重力異常數(shù)據(jù)在一定波段內(nèi)以線性關系為主，疊加非線性二次項和三次項可較好地提升反演精度，高于三次項的部分對反演的結果沒有明顯的作用。

2 衛(wèi)星測高數(shù)據(jù)反演海底地形的原理及方法

2.1 衛(wèi)星測高數(shù)據(jù)計算海洋重力異常

基于Sandwell方法，可從衛(wèi)星測高數(shù)據(jù)得到沿軌海面高程斜率，結合測高點的位置和時間信息得到沿軌大地水準面起伏對空間的導數(shù)，從而計算海面垂線偏差。如圖1所示，以相對于平均海面z0深度的參考水深面作xy平面建立直角坐標系，x軸指向地理東向，y軸指向地地理北向，z軸與x和y軸構成右手坐標系。記海面東向和北向垂線偏差的分量分別為

(1)

(2)

式中：N(x,y)表示點(x,y)處大地水準面高度，即大地水準面與參考橢球面的距離差。

圖1 平均海面與海底地形坐標示意圖Fig.1 Coordinate system for mean sea surface and sea floor topography

根據(jù)逆Vening-Meinesz公式，可由垂線偏差計算海洋重力異常Δg[11]為

η(P)sinα]dσ

(3)

式中：g0為海平面的重力加速度，通常取9.8m·s-2；ξ(P)和η(P)分別為P處垂線偏差的北向和東向分量；α為流動點至測量點的方位角；H′(ψ)為積分核函數(shù)。

(4)

式中：ψ為流動點到計算點的球面角距。

至此，已建立測高數(shù)據(jù)與海洋重力異常的相互關系。根據(jù)物理大地測量學相關理論，海底地形起伏變化將引起局部重力場擾動，記為擾動重力δg，該重力場擾動包含于海洋重力異常。針對重力異常數(shù)據(jù)在相關波段內(nèi)濾波，保留與海底地形具有高相關性的擾動重力，才可進而反演海底地形。

2.2 利用擾動重力反演海底地形

假設存在z=h(r)的一個海底地形起伏面，如圖1所示，h(r)為相對于參考水深面的海底地形高度。圖中海面上P點處海底地形產(chǎn)生的擾動位記為

(5)

式中：G=6.672 6×10-11Nt·m2·kg-2為地球引力常數(shù)；r0和r分別為計算點的位置矢量和海底地形起伏面z=h(r)的流動體元的位置矢量；Δρ為海底固體巖石層與海水的密度差。

采用頻域分析法，較之經(jīng)典的時域分析法可以快速地體現(xiàn)系統(tǒng)的響應過程，高效率地計算得到系統(tǒng)響應輸出，并為相關波段的截取提供便利。將式(5)在z=z0平面內(nèi)進行二維傅里葉變換可快速得到

(6)

海面擾動重力是重力擾動位U(r0)在垂直方向上的導數(shù)，計算式(6)并根據(jù)傅里葉變換的性質可得頻域內(nèi)擾動重力δg的表達式：

F[δg]=|k|F[U(r0)]=

(7)

根據(jù)式(7)即可利用擾動重力δg數(shù)據(jù)反演海底地形高度h(r)[1,12-13]。

2.3 反演方法

本文基于式(7)，利用規(guī)則的衛(wèi)星測高重力異常數(shù)據(jù)結合稀疏的船測水深數(shù)據(jù)反演海底地形。對測高重力異常數(shù)據(jù)在與海底地形有高相關性的波段進行濾波，得到擾動重力。文中以部分船測點作為控制點，將控制點處的水深數(shù)據(jù)作為真實水深數(shù)據(jù)，濾除與擾動重力相關性差的長波水深數(shù)據(jù)得到殘余水深。以船測殘余水深數(shù)據(jù)與船測航跡上的擾動重力建立相關關系，得到只考慮一階項的線性系數(shù)和考慮多階項的非線性系數(shù)。在船測航跡之外的區(qū)域，將得到的線性和非線性系數(shù)疊加到擾動重力的作用中，從而得到非船測航跡處預測水深。結合預測水深數(shù)據(jù)與長波水深數(shù)據(jù)，即可得到計算區(qū)域內(nèi)的全波段水深，即海底地形。將未參與計算的船測點作為外部校核點，用以衡量海底地形的反演精度。

3 數(shù)據(jù)準備

仿真實驗針對區(qū)域40°N～50°N，150°E～160°E的北太平洋海域內(nèi)，基于一組衛(wèi)星測高重力異常數(shù)據(jù)進行海底地形反演，并應用一組稀疏船測水深數(shù)據(jù)作為控制點輔助仿真。重力異常數(shù)據(jù)來源于丹麥科技大學空間研究所，數(shù)據(jù)分辨為1′×1′，數(shù)據(jù)總點數(shù)為361 201，該數(shù)據(jù)由多年的測高衛(wèi)星數(shù)據(jù)計算得到，數(shù)據(jù)如圖2所示。

圖2 重力異常數(shù)據(jù)Fig.2 Gravity anomalies data

稀疏船測數(shù)據(jù)來源于美國大氣海洋管理局(NOAA)，船測數(shù)據(jù)總點數(shù)為57 986。稀疏船測水深數(shù)據(jù)(如圖3所示)表示測量船(海面)到海底地形起伏面的深度，取向上為正。

圖3 船測水深數(shù)據(jù)(圖中包含測量船的航跡線)Fig.3 Shipboard bathymetry data(with survey ship routes)

文中取56 826個船測數(shù)據(jù)點作為控制點帶入海底地形反演的計算，數(shù)據(jù)點分布如圖4所示。

圖4 船測水深控制點數(shù)據(jù)(無船測航跡線)Fig.4 Shipboard bathymetry data of reference points(without survey ship routes)

其余1160個船測數(shù)據(jù)點作為核驗點帶入海底地形反演精度的核驗，數(shù)據(jù)點分布如圖5所示。

圖5 船測水深核驗點數(shù)據(jù)Fig.5 Shipboard bathymetry data of check points

針對船測水深數(shù)據(jù)中部分數(shù)據(jù)偏差較大的情況，需要對船測水深數(shù)據(jù)進行粗差校正。為了方便粗差校正，并濾除長波水深數(shù)據(jù)，首先對不規(guī)則的船測水深數(shù)據(jù)插值得到格網(wǎng)化的原始船測水深，如圖6所示；經(jīng)過粗差校正后的水深數(shù)據(jù)如圖7所示，從圖6與圖7的對比可看出，粗差校正有效地剔除船測水深中偏差較大的數(shù)據(jù)，從而提高仿真實驗的數(shù)據(jù)可信性。

圖6 船測水深格網(wǎng)化后的三維圖像Fig.6 3D presentation for gridded shipboard bathymetry

圖7 粗差校正后的船測水深三維圖像Fig.7 3D presentation for gridded shipboard bathymetry after cross error correction

4 仿真分析

仿真基于文中第2節(jié)的反演原理和方法，利用第3節(jié)的實驗準備數(shù)據(jù)進行仿真實驗，具體的仿真流程如圖8所示。

圖8 仿真實驗流程圖Fig.8 Flow chart for simulation

根據(jù)重力異常與水深在海底地形波長為20～200 km波段內(nèi)具有高相關性的理論，對船測水深數(shù)據(jù)進行低通濾波，濾波器設計如圖9所示，可獲得波長大于200 km的長波船測水深數(shù)據(jù)，如圖10所示(圖中將船測水深轉化為地形高度)；同時對重力異常數(shù)據(jù)進行帶通濾波，濾波器設計如圖11所示，可獲得波長介于20 km與200 km之間的擾動重力數(shù)據(jù)并向下延拓至參考水深處(文中參考水深為z0=4 km)，如圖12所示。

圖9 低通濾波器Fig.9 Low pass filter

圖10 長波船測海底地形數(shù)據(jù)Fig.10 Long wavelength shipboard bathymetry

圖11 帶通濾波器Fig.11 Band pass filter

針對船測殘余水深數(shù)據(jù)與20～200 km的船測航跡上的擾動重力數(shù)據(jù)進行回歸分析并均勻插值，得到船測航跡上海底地形與擾動重力的線性系數(shù)和多階非線性系數(shù)。將仿真區(qū)域的擾動重力與相關系數(shù)相乘即可得到短波水深數(shù)據(jù)，將計算水深與長波水深相加即可得到全區(qū)域的反演水深。

圖12 帶通濾波并向下延拓至參考水深處的擾動重力數(shù)據(jù)Fig.12 Gravity anomaly disturbance data being band pass filtered and downward continued to reference water depth

經(jīng)過計算，最終得到1′×1′的海底地形模型。圖13所示為僅考慮一階線性項的海底地形反演結果，地形高度以相對于參考水深面(低于平均海面4 km處的參考面)的地形起伏高度進行表示。由于考慮非線性項與只考慮線性項的反演結果差異遠小于海底地形的量級，因此考慮非線性項與只考慮線性項的反演結果圖沒有明顯區(qū)別，文中只給出取線性項的海底地形反演圖，具體非線性項的影響在表1中以數(shù)據(jù)形式體現(xiàn)。

圖13 取線性項的海底地形反演結果Fig.13 Predicted bathymetry with first order

反演的精度可由反演水深結果與未參與計算的船測控制點水深進行比對，表1對考慮一階線性項和多階非線性項的反演數(shù)據(jù)精度進行統(tǒng)計。

表1 海底地形反演數(shù)據(jù)統(tǒng)計

將反演得到的海底地形與船測數(shù)據(jù)進行相關性檢驗可得數(shù)據(jù)如圖14所示。

圖14 反演海底地形與船測數(shù)據(jù)的相關性Fig.14 Coherence between predicted bathymetry and shipboard bathymetry

根據(jù)表1及圖14可看出，在20～200 km的波段范圍內(nèi)，重力異常與海底地形以線性關系為主，疊加非線性二次項和三次項將反演精度的平均值量級提升至米級，反演得到的海底地形與船測數(shù)據(jù)的相關性也有較大提高，高于三次項的部分對反演的結果沒有明顯的作用，同時由于控制點數(shù)據(jù)量有限，更高階的關系系數(shù)并不能準確地獲得，從而導致從五階之后階數(shù)增大反演結果卻并不理想。

5 結束語

本文詳細研究并推導了由衛(wèi)星測高數(shù)據(jù)反演海底地形的原理及方法，基于此方法，利用衛(wèi)星測高重力異常數(shù)據(jù)和船測水深數(shù)據(jù)對北太平洋10°×10°的區(qū)域進行海底地形反演。數(shù)據(jù)及分析結果表明：該方法可以有效地反演海底地形，且海底地形與重力異常數(shù)據(jù)在20～200 km波段內(nèi)主要呈現(xiàn)線性關系，非線性二次項和三次項對反演精度的提升具有一定的作用，高于三次項的部分對反演的結果沒有明顯的作用。由于波長小于20 km和波長大于200 km海底地形起伏引起重力異常的物理機理非常復雜，且重力異常在向下延拓時會出現(xiàn)解不唯一且不穩(wěn)定的問題，由此產(chǎn)生的低于20 km和高于200 km波段的重力異常數(shù)據(jù)不能準確地與船測水深數(shù)據(jù)建立關系，也是導致最終海底地形反演精度不高的原因之一。因此，未來對全波段重力異常與水深關系的研究將具有非常重要的意義。

References)

[1]Walter H F Smith, David T Sandwell. Bathymetry prediction from dense satellite altimetry and sparse shipboard bathymetry[J]. Journal of Geophysical Research, 1994, 99(B11): 21803-21824

[2]T H Bell. Statistical features of sea-floor topography[J]. Deep Sea Research, 1975, 22: 883-892

[3]D C McAdoo, K M Marks. Gravity fields from Geosat and ERS-1 satellite altimetry for geodynamic applications[J]. Eos Trans, 1992(74): 133

[4]D T Sandwell, Walter H F Smith. Global marine gravity from ERS-1, Geosat and Seasat reveals new tectonic fabric[J]. Eos Trans,1992, AGU(73):133

[5]G A Neumann, D W Forsyth, D T Sandwell. Comparison of marine gravity from shipboard and high-density satellite altimetry along the mid-Atlantic ridge, 30.5°-35.5°S[J]. Geophysical Research Letters, 1993, 20(15): 1639-1642

[6]呂華慶. 物理海洋學基礎[M]. 北京：海洋出版社, 2012

LU Huaqing. Fundamental of physical oceanography[M]. Beijing: Ocean Press, 2012 (in Chinese)

[7]T H Dixon. Bathymetric predition from Seasat altimeter data[J]. Journal of Geophysics and Research, 1983, 88: 1563-1571

[8]G Ramillien, A cazenave. global bathymetry derived from altimeter data of the ERS-1 geodetic mission[J]. Journal of Geodynamics, 1997, 23: 129-149

[9]C Hwang. A bathymetric model for the South China Sea from satellite altimetry and depth data[J]. Marine Geodesy, 1999, 22(1): 37-51

[10]黃謨濤,翟國君，歐陽永忠, 等. 衛(wèi)星測高資料在反演海底地形中的應用[J]. 海洋測繪, 2002, 22(1): 3-7

Huang Motao, Zhai Guojun, Ouyang Yongzhong, et al. Applications of satellite altimetry data in bathymetry prediction[J]. Hydrographic Surveying and Charting, 2002, 22(1): 3-7 (in Chinese)

[11]C Hwang. Inverse Vening-Meinesz formula and deflection-geoid formula: applications to the predictions of gravity and geoid over the South China Sea[J]. Journal of Geodesy, 1998, 72(5): 304-312

[12]Hu Minzhang, Li Jiancheng, Li Hui, et al. Bathymetry predicted from vertical gravity gradient anomalies and ship soundings[J]. Geodesy and Geodynamics, 2014, 5(1): 41-46

[13]歐陽明達,孫中苗,翟振和, 等.采用重力異常的導納理論推估海底地形[J]. 測繪學報, 2015, 44(10): 1092-1099

Ouyang Mingda, Sun Zhongmiao, Zhai Zhenhe, et al. Bathymetry prediction based on the admittance theory of gravity anomalies[J]. Acta Geodaetica et Cartographica Sinica, 2015, 44(10): 1092-1099 (in Chinese)

(編輯：張小琳)

Bathymetry Prediction and Analysis Based on Satellite Altimetry Data

SUI Xiaohong1,2ZHANG Running3WAN Xiaoyun2LI Yang1

(1 DFH Satellite Co., Ltd., Beijing 100094, China) (2 Qian Xuesen Laboratory of Space Technology, Beijing 100094, China) (3 Beijing Institute of Spacecraft System Engineering, Beijing 100094, China)

Both of the resolution and accuracy of bathymetry model is far from the model of the earth land and geoid, because of observing difficulty. The traditional measurement of the bathymetry is using shipboard sounding technique, which has the limitations of low efficiency, sparse and nonuniform covering. In this paper, bathymetry prediction method based on satellite altimetry data is presented. An 1′×1′ bathymetry model of a region in north Pacific is predicted through simulation in this paper, and the nonlinear relationship between satellite altimetry data and bathymetry data is analyzed. The results show that, it is feasible to use this method predicting bathymetry. The relationship between satellite altimetry gravity anomalies and bathymetry is mainly linear in certain wavelength; adding the second and third order nonlinear terms will improve the accuracy of the inversion, however the higher order terms have little effect.

bathymetry; satellite altimetry; gravity anomalies; shipboard sounding technique

2017-02-01；

2017-05-26

國家自然科學基金(41404019，41674026)，中國科學院太空應用重點實驗室開放基金(CSU-WX-A-KJ-2016-044)

眭曉虹，女，博士研究生，研究方向為衛(wèi)星測高總體設計及應用。Email：yuanyuansui@hotmail.com。

P228

A

10.3969/j.issn.1673-8748.2017.03.020