多种重力场数据混叠的中国海岸带海域重力似大地水准面精化

1、多种重力场数据混叠的中国海岸带海域重力似大地水准面精化国家自然科学基金（41374081），国家科技支撑计划课题（2012BAB16B01）第一作者: 郭春喜，男，2008年获得武汉大学大地测量学与工程测量博士学位，主要从事大地测量数据处理与地球重力场研究；通讯作者：章传银，男，出生于1968年12月，主要从事大地测量与地球重力场研究，zhangchy郭春喜1,章传银2,王斌1,柯宝贵2,王夏莉1,蒋涛2,范宏涛1,王伟21国家测绘地理信息局大地测量数据中心，西安 7100542中国测绘科学研究院，北京100830摘要：本文采用统一的地球外部（含地面）任意类型重力场参数的精密归算方案，对中国沿

2、海地区陆地重力、海岸带航空重力、海域船测重力及卫星测高等多种重力场数据进行了精细化处理，进而以大地水准面为边界面，按不同来源重力数据的可靠性进行配权，采用加权最小二乘谱分析方法实现了陆海重力场的数据集成，求得中国海岸带海域陆海统一的2.52.5地面扰动重力数值模型，最后，在Stokes框架下，采用严密的球近似算法，计算了中国海岸带海域陆海统一的2.52.5重力似大地水准面数值模型。多种方案的检测分析结果显示，2.52.5地面扰动重力精度达到毫伽级，重力似大地水准面精度达到厘米级。关键词：多种数据混叠；海岸带；重力似大地水准面；Stokes框架Refinement Gravimetric Qua

3、si-geoid Based on Various Gravity Data Hybrid in Chinese Coastal ZoneGUO Chun-Xi1,ZHANG Chuan-Yin2,WANG Bin1,KE Bao-Ggui2,WANG Xia-Li1,JIANG Tao2,FAN Hong-Tao1 ,WANG Wei2 ,MA Xin-Ying1,LI Wei21 Geodetic Data Processing Center, State Bureau of Surveying and Mapping, Xian 710054, China2 Chinese Academ

4、y of Surveying and Mapping, Beijing 100830, ChinaAbstract: In this paper, the off-shore terrestrial gravity data, coastal zone airbornegravity data , Marine ship gravity data, satellite altimetry gravity data and so on ,have been processed based on an unified accurate reduction scheme, which weresui

5、ted to any gravity field parameters outside or including the ground of the earth. the various gravity data took geoid as a boundary surface and were assigned to different weight according to their quality reliability. An weighted least square spectrum analysis method was adopted to implementedland a

6、nd ocean gravity data integration. Then a gravity disturbance numerical model was obtained, with 2.5 minutes spatial resolution, in region of chinese coastal zone. Under the stokes frame the gravimetric Marine Quasi-geoid, which spatial resolution was also 2.5mimutes, was calculated through the rigo

7、rous spherical approximation. Several testing schemes were compared and the results shown that: gravity disturbance numerical model and gravimetric Quasi-geoid can be achieved mgal and centimeter level respectively.Key Word: Various Gravity Data Hybrid; Coastal Zone; gravimetric Quasi-geoid; Stokes

8、Frame1引言本文的重力场及似大地水准面精化区域覆盖中国海岸带、海域及沿海省市，精化区域内重力数据源复杂（图1），主要有：陆地地面重力、海面船测重力、飞行高度上的航空重力和海域卫星测高等数据，这些重力资料的典型特征是：数据分布不规则，类型多样，可靠性（精度）差异大，数据获取高度不同，存在重力空白区。这种异常复杂情况下的局部重力场逼近及大地水准面精化是当前国际上的热点和难点问题12。考虑到海岸带海域地区的地面海拔高普遍较小，因此我们选择大地水准面作为统一的边界面，对离散的地面重力、船测重力、航空重力和浅水卫星测高采用了完全一致的地形影响和重力归算方法，以提高多源重力场数据融合和局部重力场逼近性

9、能3。重力似大地水准面计算采用“局部地形影响+参考重力场”组合移去恢复法，参考重力场模型选用融合了Lageos、GRACE、地面重力、卫星测高和最新GOCE数据的EIGEN6C2（1949阶）重力场模型，截断到720阶，积分半径30。移去地形影响时，采用1515的数字高程模型DEM（海域高程置零），恢复地形影响时采用2.52.5数字高程模型，地形影响计算的积分半径为60。为满足参考重力场移去恢复法的适用条件，减少重力似大地水准面计算过程的积分运算次数， 本文在Stokes框架中，由大地水准面上的扰动重力按严密的球近似算法直接一步计算地面重力高程异常，得到重力似大地水准面数值模型4。为便于计算，

10、我们将精化区域分成如下四个区域：山东测区、浙江测区、海南测区和南海测区，如图1。测区之间重合2，拼接时要求2.52.5地面扰动重力互差的标准差小于0.1mGal，重力似大地水准面互差的标准差小于0.1cm。计算软件采用中国测绘科学研究院研发的高精度局部重力场计算平台PALGrav2.0。2离散重力数据处理本文的离散重力场数据处理目的是，将实测重力数据归算到大地水准面上，且要求归算后的重力场参数，不论是地面重力、船测重力、航空重力还是测高重力场，也不论实测点是在陆地还是在海域，都应是同一物理量，以便下一步进行陆海重力场集成。2.1统一的数据处理算法（1）计算实测点的离散扰动重力当实测点的高程为正

11、常高时，则需利用EIGEN6C2模型和高程基准差别（章传银等，2009），将正常高转换为大地高；海面大地高则由高精度测高平均海面高模型（MSS）内插得到。由实测重力和大地高计算测点的扰动重力，得到地面点扰动重力、航空测线扰动重力、船载测线扰动重力。对于卫星测高数据，先反演测高扰动重力，再离散化后求得。（2）全部测点进行局部地形影响改正。局部地形影响采用严密的球近似算法。航空高度上的局部地形影响也直接计算（章传银等，2009）；海面测点（船载测线扰动重力和测高扰动重力）的局部地形影响与陆地地面局部地形影响算法完全相同。本文不进行海水地形改正。（3）当实测点正常高不为零时，进行解析延拓改正地面点和

12、航空测线高度上的扰动重力需要解析延拓到大地水准面（海面）。为提高延拓的稳定性，延拓前先移去局部地形影响；为大幅抑制边缘效应，本文还移去了EIGEN6C2（1949阶）满阶的模型扰动重力，下称模型延拓改正；最后对扣除局部地形影响和模型扰动重力的残差扰动重力按1阶梯度法（积分半径最大10）完成由测点高度到大地水准面（海面）的解析延拓，下称剩余延拓改正。（4）当实测点或其正下方的地形高（正常高）不为零时，进行布格改正地面实测扰动重力点和航空重力测线落到陆地的测线扰动重力需要进行布格改正，布格改正的地形高采用1515DEM内插值。测试表明，海域扰动重力代表性误差普遍小于海域布格扰动重力代表性误差，因此

13、，本文不进行海水布格改正。（5）陆地布格扰动重力与海域剩余扰动重力本文陆地和海域测点的重力改正项不同。为便于说明问题，这里暂且将陆地扰动重力经过上述改正后称为布格扰动重力，将海域扰动重力经过上述改正后称为海域剩余扰动重力。不难看出，海域剩余扰动重力与陆地布格扰动重力是大地水准面上的同一物理量，在陆海交界处无缝拼接，本文进一步将陆海统一的布格（剩余）扰动重力简称为布格（剩余）扰动重力。2.2地面离散重力数据处理地面离散重力数据处理主要步骤包括：测点扰动重力计算，局部地形影响计算，平面布格改正计算、模型延拓改正计算（EIGEN6C2满阶模型扰动重力移去恢复法），实测数据粗差探测，残差扰动重力剩余延

14、拓改正计算，大地水准面上扰动重力和布格扰动重力计算。表1地面离散重力数据处理过程统计（mGal）Table 1 statistic of the terrain scatter gravity data of different procedure (mGal)序号类型最大值最小值平均值标准差1实测扰动重力165.789-218.804-14.76723.5893局部地形影响0.086-34.424-1.0241.7323平面布格改正0-425.018-42.07457.3854模型延拓改正92.606-56.0870.1992.6465剩余延拓改正71.660-42.5390.0150.94

15、26大地水准面扰动重力217.775-232.402-14.56824.7927大地水准面布格扰动重力65.718-349.847-55.60456.5282.3船载重力测线数据处理船载重力测线数据处理主要步骤包括：测线扰动重力计算，局部地形影响计算，测线粗差探测及平差计算，海域剩余扰动重力计算。分析发现，我们收集到的中国近海船测重力数据，不同期测量和不同单位测量成果在重叠区域存在不一致现象，为有效地消除原始数据中的各种矛盾，提高船测数据处理精度，本文采用如下数据处理方案对船测数据进行了特殊处理：首先基于扰动重力数据逐一对每条测线数据进行沿线粗差探测；其次，基于测线交叉点处的扰动重力不符值，优

16、选几个高精度的参考测线；再次，基于参考测线对区域的测线网进行平差，从而消除不同区域的系统偏差；最后，将处理前后的船测扰动重力与测高扰动重力独立进行比较，评价数据处理方案的可靠性。表2纠正前后部分测区船载测线交叉点不符值变化（mGal）Table 2 Comparison of the bias value at the cross point tetween pre-rectify and post-rectify ship gravity survey line in the part region offshore of China （mGal）测区编号最大值最小值平均值标准差1纠正前5.

17、415-5.940-0.0601.842纠正后0.965-2.068-0.0220.3082纠正前5.781-5.126-0.0862.030纠正后2.190-1.856-0.0010.5133纠正前4.468-4.911-0.1161.788纠正后1.409-1.167-0.0160.2804纠正前5.086-4.976-0.2832.645纠正后2.543-2.2640.0060.7865纠正前22.583-13.1290.1223.394纠正后2.273-2.339-0.0020.2632.4航空重力测线数据处理航空重力测线数据处理主要步骤包括：测线扰动重力计算，测线高度局部地形影响计算

18、，陆地平面布格改正计算、模型延拓改正计算，测线粗差探测及平差计算，测线偏差纠正，残差扰动重力剩余延拓改正计算，大地水准面上布格（剩余）扰动重力计算。陆地航空重力测线扰动重力归算与离散地面扰动重力归算方法完全相同，航空重力测线粗差探测及平差计算与船载测线粗差探测及平差计算完全相同。考虑到航空重力数据主要分布在陆地、陆海交界空白区及海域上空，为检核或解决参考测线选择不可靠导致测线重力值总体偏移问题，且考虑到陆地重力精度较高，本文采用与航空重力重合的陆地数据对平差后的航空重力进行剩余偏差纠正；并用与航空重力测量区域重合的船测数据对航空重力数据处理结果独立进行检核。表3纠正前后部分测区船测与航空扰动重

19、力之差统计（mGal）Table 3 Statistic of the gravity disturbance discrepancy value between ship gravity and airborne gravity of pre-rectify/ post-rectify in the part region offshore of China （mGal）区号状态最大值最小值平均值标准差1纠正前11.291-11.888-0.8484.376纠正后12.291-10.8880.1524.3762纠正前79.070-46.0053.4586.980纠正后76.370-48.70

20、50.7586.9803纠正前9.435-2.2704.1912.214纠正后7.281-4.4242.0372.2144纠正前12.959-7.5591.6393.815纠正后11.599-8.9190.2793.8152.5浅水测高边界流探测与抑制在浅海水域，卫星测高质量偏低，海面地形复杂，存在多种不同尺度的边界流，卫星测高数据反演重力场的可靠性明显降低。本文采用重力场积分迭代法来探测和改善浅水卫星测高重力场的质量，取得了较好的效果。图2（单位：mGal）中，蓝色标注值为卫星测高反演的原始残差扰动重力（扰动重力局部地形影响Eigen720阶模型扰动重力），黑色标注值为陆地残差布格扰动重力（

21、布格扰动重力Eigen720阶模型扰动重力），红色标注值为改善后的卫星测高残差剩余扰动重力。图2浅水测高重力改善前后与陆地重力衔接情况对照Fig 2 Comparison the terrain gravity data and the shallow water region altimetry gravity data before or after correction由图2可以看出，改善前一些陆海交界处测高重力与陆地重力不能平滑衔接，有的相差超过30mGal，改善后陆海交界处陆地与测高数据的一致性得到明显改善。由此可以认为，经上述处理后，浅水卫星测高重力场的可靠性得到一定程度的改善，沿

22、岸边界流影响得到抑制。表4边界流抑制前后测高残差剩余扰动重力统计（mGal）Table 4 Statistics ofthe residualsdisturbance gravity of the altimetry data before or after the boundary current suppressing (mgal)测区最大值最小值平均值标准差山东测区改善前23.734-28.457-0.2015.593改善后17.406-14.654-0.1523.735浙江测区改善前31.303-40.212-0.0685.518改善后12.302-17.2680.0423.613在卫

23、星测高数据处理过程中，我们探测出三处沿岸狭长的边界流，他们分别是河北辽宁交界沿岸的渤海边界流、山东沿岸的黄海边界流和福建沿岸的东海边界流。3陆海多源重力场数据集成陆海重力场数据集成（局部重力场逼近）的主要目标是：基于多源重力场离散数据处理后的大地水准面上布格（剩余）扰动重力数据，计算2.52.5陆海一致的布格（剩余）扰动重力数值模型。由于多源重力场数据在陆海交界处分布不规则，类型多样，数据的可靠性（精度）差异较大，且存在较多的重力空白区（特别是陆海交界的浅水海域），为区别对待不同来源的重力场数据，本文采用加权最小二乘谱分析方法进行多源重力场数据集成。按数据源的可靠性由高到低排序，分别对大地水准

24、面上的陆地、船测、航空和测高布格（剩余）扰动重力据配置不同的权值，采用最小二乘谱分析方法进行格网化，得到陆海一致的2.52.5布格（剩余）扰动重力格网。由陆海一致的布格（剩余）扰动重力格网内插离散重力点的布格（剩余）扰动重力，并计算残差，分别统计陆地、船测、航空、测高等扰动重力的代表性误差。根据不同类型扰动重力代表性误差大小重新考察权比关系的合理性，调整不合理的配权，重新进行最小二乘谱分析法重力场集成，求得最终的2.52.5陆海一致布格（剩余）扰动重力数值模型。表5浙江测区陆海重力代表性误差统计（mGal）Table 5 Statistics of the land and seagravit

25、yrepresentativeerror in Zhejiang（mGal）数据类型最大值最小值平均值标准差陆地27.798-40.966-3.3804.460船测36.021-20.6520.0892.738航空17.066-18.117-1.5883.413测高44.641-35.789-0.3903.172国际船测41.943-35.3831.0025.5474地面扰动重力与重力似大地水准面计算令残差剩余扰动重力=扰动重力-扰动重力局部地形影响-模型扰动重力。本文采用完全相同的计算流程（如图3）由大地水准面扰动重力计算2.52.5地面扰动重力和重力似大地水准面数值模型。地面扰动重力计算时

26、采用Possion积分公式，地面重力高程异常（重力似大地水准面）计算时采用广义Hotine积分公式5。算法采用“局部地形影响+模型重力场”组合移去恢复法。参考重力场模型EIGEN6C2截断到720阶，地形影响恢复时采用2.52.5DEM（海面置零）。Possion积分和广义Hotine积分的积分半径采用30，局部地形影响的积分半径采用60。大地水准面上2.52.5残差剩余扰动重力布格改正参考重力场模型2.52.5布格扰动重力待计算参数局部地形影响计算重力场广义积分运算待计算参数地面模型值计算2.5数字高程模型2.52.5地面重力场数值模型待估参数地面残差剩余量图3地面重力场数值模型计算流程Fi

27、g 3 The calculating flow chart of the ground gravitynumericalmodel4计算结果检核与分析4.1平均扰动重力代表性误差估计以陆地重力、船测重力、航空重力和测高重力分类，分别按表5统计山东、浙江、海南和南海四个测区的扰动重力代表性误差，结果如表6。表6中，海域代表性误差以数据量大或代表性误差大的类型数据统计结果代替。表6陆海平均扰动重力代表性误差统计Table 6 Statistics of the land and seamean Gravity disturbancerepresentativeerror测区名称陆地mGal海域m

28、Gal备 注山东测区4.253.84海域航空重力数据量多，取航空重力代表性误差浙江测区4.462.74海域船测重力数据量多，取船测重力代表性误差海南测区6.193.11海域测高数据量多，取测高重力代表性误差南海测区无数据5.33南海国际船测代表性误差最大，取国际船测重力代表性误差表6中海南测区陆地扰动重力代表性误差最大，达到6.2mGal，是因为本文能收集的海南测区陆地重力数据稀少导致的，如图1。总体统计结果显示：本文的2.52.5平均扰动重力代表性在陆地为4.6mGal，海域4.1mGal。4.2由代表性误差推估重力似大地水准面误差重力似大地水准面误差N与平均扰动重力代表性误差g之间的关系可

29、表示为如下正比关系68：（cm）（mGal）式中：为空间分辨率大小，单位为分。利用上式，采用表6中2.52.5扰动重力代表性误差推估陆海2.52.5重力似大地水准面误差，结果如表7。表7重力代表性误差推算的重力似大地水准面误差Table 7 Gravimetric Quasi-geoid error derived from the gravity representativeerror测区名称扰动重力代表性误差mGal重力似大地水准面误差cm陆地海域陆地海域山东测区4.253.843.042.75浙江测区4.462.743.191.96海南测区6.193.114.432.22南海测区无数据5

30、.33无数据3.814.3重力似大地水准面外部检核（1）利用沿岸13处海岛GNSS三角高程跨海传递成果，对重力似大地水准面数值模型进行外部检核，计算残差GNSS水准高程异常，统计结果为：标准差3.7cm , 最大值-1.7cm，最小值-18.4cm，平均值-12.4cm。（2）利用中国国家测绘地理信息局区域似大地水准面精化项目布设的106个沿岸陆地A、B级大地控制点GPS/水准成果，对重力似大地水准面数值模型进行外部检核，计算残差GNSS水准高程异常，统计结果为：标准差3.9cm , 最大值-3.6cm，最小值-23.4cm，平均值-12.1cm。（3）利用收集到的400个沿海省市陆地GNSS水准点对重力似大地水准面数值模型进行外部检核，计算残差GNSS水准高程异常，统计结果为：标准差4.2cm，最大值1