留言板

尊敬的读者、作者、审稿人, 关于本刊的投稿、审稿、编辑和出版的任何问题, 您可以本页添加留言。我们将尽快给您答复。谢谢您的支持!

姓名
邮箱
手机号码
标题
留言内容
验证码

综合半参数核估计与正则化方法的航空重力向下延拓模型分析

伍丰丰 黄海军 任青阳 樊文有 陈洁 潘雄

伍丰丰, 黄海军, 任青阳, 樊文有, 陈洁, 潘雄. 综合半参数核估计与正则化方法的航空重力向下延拓模型分析[J]. 武汉大学学报 ● 信息科学版, 2020, 45(10): 1563-1569. doi: 10.13203/j.whugis20180491
引用本文: 伍丰丰, 黄海军, 任青阳, 樊文有, 陈洁, 潘雄. 综合半参数核估计与正则化方法的航空重力向下延拓模型分析[J]. 武汉大学学报 ● 信息科学版, 2020, 45(10): 1563-1569. doi: 10.13203/j.whugis20180491
WU Fengfeng, HUANG Haijun, REN Qingyang, FAN Wenyou, CHEN Jie, PAN Xiong. Analysis of Downward Continuation Model of Airborne Gravity Based on Comprehensive Semi-parametric Kernel Estimation and Regularization Method[J]. Geomatics and Information Science of Wuhan University, 2020, 45(10): 1563-1569. doi: 10.13203/j.whugis20180491
Citation: WU Fengfeng, HUANG Haijun, REN Qingyang, FAN Wenyou, CHEN Jie, PAN Xiong. Analysis of Downward Continuation Model of Airborne Gravity Based on Comprehensive Semi-parametric Kernel Estimation and Regularization Method[J]. Geomatics and Information Science of Wuhan University, 2020, 45(10): 1563-1569. doi: 10.13203/j.whugis20180491

综合半参数核估计与正则化方法的航空重力向下延拓模型分析

doi: 10.13203/j.whugis20180491
基金项目: 

国家自然科学基金 41631072

国家自然科学基金 41874009

国家自然科学基金 11471105

详细信息
    作者简介:

    伍丰丰,硕士生,主要研究方向为航空重力数据处理。WuFF@cug.edu.cn

    通讯作者: 潘雄,博士,教授。pxjlh@163.com
  • 中图分类号: P223

Analysis of Downward Continuation Model of Airborne Gravity Based on Comprehensive Semi-parametric Kernel Estimation and Regularization Method

Funds: 

The National Natural Science Foundation of China 41631072

The National Natural Science Foundation of China 41874009

The National Natural Science Foundation of China 11471105

More Information
    Author Bio:

    WU Fengfeng, postgraduate, specializes in airborne gravity data processing. E-mail: WuFF@cug.edu.cn

    Corresponding author: PAN Xiong. PhD, professor. E-mail:pxjlh@163.com
  • 摘要: 在航空重力向下延拓过程中,将重力数据中的系统误差和离散化造成的模型误差用非参数分量表达。在无外部数据的情况下,建立基于半参数核估计方法的重力向下延拓模型,为了改善泊松积分离散后的设计矩阵的病态影响,引入正则化方法,提出了综合半参数核估计和正则化方法的逆泊松积分延拓方法。基于EGM2008(earth gravity model 2008)模型计算了某地空中重力异常,采用线性项和周期项系统误差进行仿真实验,以及美国某地实测重力异常数据,验证了本文方法在改善病态性和分离系统误差方面的有效性。结果表明,本文方法在无外部数据时,能有效地分离系统误差并具有较高的精度。
  • 图  1  向下延拓值与基准值的差值统计

    Figure  1.  Differences Between the Results of Downward Continuation and Real Data

    图  2  两方案6.3 km向下延拓系统误差估计情况

    Figure  2.  Systematic Errors Estimation of 6.3 km Downward Continuation of Two Schemes

    表  1  仿真区域各项重力异常特性的统计

    Table  1.   Characteristics of Gravity Anomalies at Simulation Area

    高度/km 最小值
    /mGal
    最大值
    /mGal
    平均值
    /mGal
    标准差
    /mGal
    0 -104.35 35.67 -29.78 32.93
    3.3 -90.93 28.03 -28.93 30.07
    6.3 -82.97 22.73 -28.19 27.84
    下载: 导出CSV

    表  2  不同方法延拓值与基准值差值

    Table  2.   Differences Between the Results and Real Data from Different Models

    延拓高度
    /km
    系统误差 随机噪声
    /km
    LS/mGal TIK/mGal SEM/mGal
    平均值 均方根 平均值 均方根 平均值 均方根
    3.3 Simu1 ±1 3.38 5.01 3.06 3.96 0.06 2.23
    ±3 3.49 10.77 2.52 4.66 0.08 3.92
    6.3 Simu1 ±1 3.49 9.65 3.20 4.41 0.06 2.78
    ±3 3.60 26.28 2.83 5.15 0.08 4.26
    3.3 Simu2 ±1 2.42 5.75 1.94 4.46 0.10 3.25
    ±3 2.47 11.02 1.70 4.95 0.10 4.35
    6.3 Simu2 ±1 2.50 9.97 2.13 5.15 0.11 3.76
    ±3 2.55 25.98 1.87 5.53 0.11 4.79
    注:Simu1、Simu2分别为方案1、方案2中的模拟系统误差
    下载: 导出CSV

    表  3  利用SEM估计观测数据中的系统误差情况

    Table  3.   Statistics of Systematic Errors Estimated by SEM Data

    延拓高度/km 系统
    误差
    随机噪声
    /mGal
    最小值
    /mGal
    最大值
    /mGal
    平均值
    /mGal
    标准差
    /mGal
    Simu1 - 1.01 5.68 3.27 1.47
    3.3 ±1 1.15 5.48 3.22 1.35
    ±3 1.37 5.51 3.31 1.32
    6.3 ±1 1.23 5.29 3.22 1.22
    ±3 1.30 5.47 3.32 1.27
    Simu2 - -3.38 7.88 2.37 3.61
    3.3 ±1 -3.57 8.04 2.32 3.14
    ±3 -3.53 7.82 2.49 3.24
    6.3 ±1 -3.38 7.98 2.32 3.14
    ±3 -3.32 7.75 2.50 3.23
    注:Simu1、Simu2分别为方案1、方案2中的模拟系统误差
    下载: 导出CSV

    表  4  无外部数据延拓值与基准值差值/mGal

    Table  4.   Differences Between the Continuation Results and Real Data from Different Models (No External Data)/mGal

    模型 最小值 最大值 平均值 均方根
    LS -29.39 33.40 2.01 6.74
    TIK -27.69 31.54 -1.42 6.43
    SEM -29.67 28.81 -0.05 6.10
    下载: 导出CSV

    表  5  平差前后交叉点不符值统计/mGal

    Table  5.   Dis crepancy Statistics of Cross Points Before and After Adjustment/mGal

    平差前后 最小值 最大值 平均值 标准差 均方根
    平差前 -6.43 5.56 -0.33 2.09 2.12
    平差后 -4.46 5.11 -0.002 1.47 1.47
    下载: 导出CSV

    表  6  部分主测线与EGM2008模型偏差统计/mGal

    Table  6.   Deviation Statistics of Part Main Survey Lines and EGM2008/mGal

    测线 偏差 测线 偏差 测线 偏差
    101 3.11 108 1.54 113 2.39
    102 2.19 109 2.46 114 1.21
    104 -0.02 110 2.46 115 2.10
    105 2.81 111 3.13 116 2.84
    106 1.76 112 1.17 117 2.53
    下载: 导出CSV

    表  7  系统误差改正后延拓值与基准值差值/mGal

    Table  7.   Difference Between the Continuation Results and the Reference Values After the Systematic Error Correction/mGal

    延拓数据 最小值 最大值 平均值 均方根
    无外部数据
    SEM改正
    -29.71 30.34 -0.05 6.10
    交叉点平
    差改正
    -27.23 30.34 1.79 6.37
    整条测试改正 -29.92 29.31 -0.26 6.12
    有外部数据
    SEM改正
    -29.34 28.38 -0.05 6.08
    下载: 导出CSV
  • [1] 孙中苗.航空重力测量理论、方法及应用研究[D].郑州: 信息工程大学, 2004

    Sun Zhongmiao. Theory, Methods and Applications of Airborne Gravimetry[D]. Zhengzhou: Information Engineering University, 2004
    [2] Bolkas D, Fotopoulos G, Braun A. On the Impact of Airborne Gravity Data to Fused Gravity Field Models[J]. Journal of Geodesy, 2016, 90(6): 561-571 doi:  10.1007/s00190-016-0893-x
    [3] Moritz H. Advanced Physical Geodesy[M]. Karlsruhe: Wichmann, 1980
    [4] 王兴涛, 石磐, 朱非洲.航空重力测量数据向下延拓的正则化算法及其谱分解[J].测绘学报, 2004, 33(1):33-38 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=chxb200401006

    Wang Xingtao, Shi Pan, Zhu Feizhou. Regularization Methods and Spectral Decomposition for the Downward Continuation of Airborne Gravity Data[J]. Acta Geodaetica et Cartographica Sinica, 2004, 33(1):33-38 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=chxb200401006
    [5] 吴太旗, 邓凯亮, 黄谟涛, 等.一种改进的不适定问题奇异值分解法[J].武汉大学学报·信息科学版, 2011, 36(8):900-903 http://ch.whu.edu.cn/article/id/636

    Wu Taiqi, Deng Kailing, Huang Motao, et al. An Improved Singular Values Decomposition Method for Ill-posed Problem[J]. Geomatics and Information Science of Wuhan University, 2011, 36(8):900-903 http://ch.whu.edu.cn/article/id/636
    [6] 蒋涛, 党亚民, 章传银, 等.基于矩谐分析的航空重力向下延拓[J].测绘学报, 2013, 42(4):475-480 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=chxb201304001

    Jiang Tao, Dang Yamin, Zhang Chuanyin, et al. Downward Continuation of Airborne Gravity Data Based on Rectangular Harmonic Analysis[J]. Acta Geodaetica et Cartographica Sinica, 2013, 42(4):475-480 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=chxb201304001
    [7] Hwang C, Shih H C, Hsiao Y S, et al. Geodetic and Geophysical Results from a Taiwan Airborne Gravity Survey: Data Reduction and Accuracy Assessment[J].Journal of Geophysical Research: Solid Earth, 2007, 112 (B4):1 548-1 554 doi:  10.1029/2005jb004220/full
    [8] Hwang C, Shih H C, Hsiao Y S, et al. Airborne Gravity Surveys over Taiwan Island and Strait, Kuroshio Current and South China Sea: Comparison of GPS and Gravity Accuracies at Different Flight Altitudes[J]. Marine Geodesy, 2012, 35(3):287-305 doi:  10.1080/01490419.2011.634962
    [9] Zhao Q L, Gabriel S, Li J C, et al. Evaluation and Comparison of the Processing Methods of Airborne Gravimetry Concerning the Errors Effects on Downward Continuation Results: Case Studies in Louisiana (USA) and the Tibetan Plateau (China)[J]. Sensors, 2017, 17(6): 1 205 doi:  10.3390/s17061205
    [10] 黄谟涛, 宁津生, 欧阳永忠, 等.联合使用位模型和地形信息的陆区航空重力向下延拓方法[J].测绘学报, 2015, 44(4):355-362 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=chxb201504001

    Huang Motao, Ning Jinsheng, Ouyang Yongzhong, et al. Downward Continuation of Airborne Gravimetry on Land Using Geopotential Model and Terrain Information[J]. Acta Geodaetica et Cartographica Sinica, 2015, 44(4): 355-362 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=chxb201504001
    [11] 石磐, 王兴涛.利用航空重力测量和DEM确定地面重力场[J].测绘学报, 1997, 26(2):117-121 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=chxb199702004

    Shi Pan, Wang Xingtao. Determination of the Terrain Surface Gravity Field Using Airborne Gravimetry and DEM[J].Acta Geodaetica et Cartographica Sinica, 1997, 26(2):117-121 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=chxb199702004
    [12] 黄谟涛, 刘敏, 邓凯亮, 等.利用重复测线校正海空重力仪格值及试验验证[J].地球物理学报, 2018, 61(8):3 160-3 169 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=dqwlxb201808006

    Huang Motao, Liu Min, Deng Kailiang.et al. Test and Correction of Scale Values for Air-sea Gravimeters Using Repeated Survey Lines[J]. Chinese Journal of Geophysics, 2018, 61 (8):3 160-3 169 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=dqwlxb201808006
    [13] Forsberg R, Olesen A V, Alshamsi A, et al. Airborne Gravimetry Survey for the Marine Area of the United Arab Emirates[J]. Marine Geodesy, 2012, 35:221-232 doi:  10.1080/01490419.2012.672874
    [14] Alberts B, Klees R. A Comparison of Methods for the Inversion of Airborne Gravity Data[J]. Journal of Geodesy, 2004, 78(1/2): 55-65 doi:  10.1002/cjg2.599/full
    [15] 欧阳永忠.海空重力测量数据处理关键技术研究[D].武汉: 武汉大学, 2013

    Ouyang Yongzhong. On Key Technologies of Data Processing for Air-sea Gravity Surveys[D].Wuhan: Wuhan University, 2013
    [16] 潘雄, 付宗堂, 范晓明.一种处理系统误差的新方法[J].仪器仪表学报, 2007 (4):630-633 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=yqyb200704010

    Pan Xiong, Fu Zongtang, Fan aoming. New Method for Processing Systematic Errors[J]. Chinese Journal of Scientific Instrument, 2007 (4):630-633 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=yqyb200704010
  • [1] 曾小牛, 李夕海, 刘继昊, 牛超, 侯维君.  一种基于改进凸集投影原理的航空重力数据插值与去噪方法 . 武汉大学学报 ● 信息科学版, 2020, 45(10): 1555-1562. doi: 10.13203/j.whugis20180470
    [2] 徐新强, 赵俊.  航空重力向下延拓的多参数正则化方法 . 武汉大学学报 ● 信息科学版, 2020, 45(7): 956-963, 973. doi: 10.13203/j.whugis20180335
    [3] 曲国庆, 孙振, 苏晓庆, 杜存鹏.  非线性参数估计的自适应松弛正则化算法 . 武汉大学学报 ● 信息科学版, 2019, 44(10): 1491-1497. doi: 10.13203/j.whugis20170425
    [4] 崔生成, 朱文越, 杨世植, 李学彬.  基于正则化的地表反射特性参数反演方法 . 武汉大学学报 ● 信息科学版, 2018, 43(8): 1264-1270. doi: 10.13203/j.whugis20160244
    [5] 陈欣, 翟国君, 暴景阳, 欧阳永忠, 陆秀平, 邓凯亮.  航空重力向下延拓的最小二乘配置Tikhonov正则化法 . 武汉大学学报 ● 信息科学版, 2018, 43(4): 578-585. doi: 10.13203/j.whugis20150728
    [6] 孙文, 吴晓平, 王庆宾, 刘晓刚, 朱志大.  基于方差分量估计的正则化配置法及其在多源重力数据融合中的应用 . 武汉大学学报 ● 信息科学版, 2016, 41(8): 1087-1092. doi: 10.13203/j.whugis20140159
    [7] 黄谟涛, 欧阳永忠, 刘敏, 翟国君, 邓凯亮.  融合海域多源重力数据的正则化点质量方法 . 武汉大学学报 ● 信息科学版, 2015, 40(2): 170-175.
    [8] 曾小牛, 刘代志, 李夕海, 苏娟, 陈鼎新, 齐玮.  一种改进的病态问题奇异值修正法 . 武汉大学学报 ● 信息科学版, 2015, 40(10): 1349-1353. doi: 10.13203/j.whugis20130709
    [9] 黄谟涛, 欧阳永忠, 刘敏, 翟国君, 邓凯亮.  海域航空重力测量数据向下延拓的实用方法 . 武汉大学学报 ● 信息科学版, 2014, 39(10): 1147-1152.
    [10] 唐利民, 朱建军.  软土路基沉降泊松模型的正则化牛顿迭代法 . 武汉大学学报 ● 信息科学版, 2013, 38(1): 69-73.
    [11] 顾勇为, 归庆明, 韩松辉, 王靳辉.  航空重力向下延拓分组修正的正则化解法 . 武汉大学学报 ● 信息科学版, 2013, 38(6): 720-724.
    [12] 高建, 谢伟, 涂志刚, 秦前清.  使用泊松方程插值方法进行遥感影像融合 . 武汉大学学报 ● 信息科学版, 2012, 37(12): 1448-1451.
    [13] 陶肖静, 朱建军, 田玉淼.  半参数模型中影响正则化参数的因素分析 . 武汉大学学报 ● 信息科学版, 2012, 37(3): 298-301.
    [14] 赵红蕊, 唐中实, 李小文.  线性正则化遥感反演中正则化参数的确定方法 . 武汉大学学报 ● 信息科学版, 2007, 32(6): 531-535.
    [15] 孙中苗, 夏哲仁, 王兴涛.  利用航空重力确定局部大地水准面的精度分析 . 武汉大学学报 ● 信息科学版, 2007, 32(8): 692-695.
    [16] 申文斌, 鄢建国, 晁定波.  重力场的局部虚拟向下延拓以及利用EGM96模型的模拟实验检验 . 武汉大学学报 ● 信息科学版, 2006, 31(7): 589-593.
    [17] 顾勇为, 归庆明, 边少锋, 郭建峰.  地球物理反问题中两种正则化方法的比较 . 武汉大学学报 ● 信息科学版, 2005, 30(3): 238-241.
    [18] 徐天河, 杨元喜.  均方误差意义下正则化解优于最小二乘解的条件 . 武汉大学学报 ● 信息科学版, 2004, 29(3): 223-226.
    [19] 徐天河, 杨元喜.  抗差Tikhonov正则化方法及其应用 . 武汉大学学报 ● 信息科学版, 2003, 28(6): 719-722.
    [20] 邓凯亮, 黄谟涛, 吴太旗, 王伟平, 欧阳永忠, 陈欣, 熊雄, 刘敏, 王许.  利用高阶径向导数带限模型进行重力向下延拓计算 . 武汉大学学报 ● 信息科学版, 0, 0(0): 0-0. doi: 10.13203/j.whugis20210630
  • 加载中
图(2) / 表(7)
计量
  • 文章访问数:  394
  • HTML全文浏览量:  93
  • PDF下载量:  28
  • 被引次数: 0
出版历程
  • 收稿日期:  2020-01-04
  • 刊出日期:  2020-10-05

综合半参数核估计与正则化方法的航空重力向下延拓模型分析

doi: 10.13203/j.whugis20180491
    基金项目:

    国家自然科学基金 41631072

    国家自然科学基金 41874009

    国家自然科学基金 11471105

    作者简介:

    伍丰丰,硕士生,主要研究方向为航空重力数据处理。WuFF@cug.edu.cn

    通讯作者: 潘雄,博士,教授。pxjlh@163.com
  • 中图分类号: P223

摘要: 在航空重力向下延拓过程中,将重力数据中的系统误差和离散化造成的模型误差用非参数分量表达。在无外部数据的情况下,建立基于半参数核估计方法的重力向下延拓模型,为了改善泊松积分离散后的设计矩阵的病态影响,引入正则化方法,提出了综合半参数核估计和正则化方法的逆泊松积分延拓方法。基于EGM2008(earth gravity model 2008)模型计算了某地空中重力异常,采用线性项和周期项系统误差进行仿真实验,以及美国某地实测重力异常数据,验证了本文方法在改善病态性和分离系统误差方面的有效性。结果表明,本文方法在无外部数据时,能有效地分离系统误差并具有较高的精度。

English Abstract

伍丰丰, 黄海军, 任青阳, 樊文有, 陈洁, 潘雄. 综合半参数核估计与正则化方法的航空重力向下延拓模型分析[J]. 武汉大学学报 ● 信息科学版, 2020, 45(10): 1563-1569. doi: 10.13203/j.whugis20180491
引用本文: 伍丰丰, 黄海军, 任青阳, 樊文有, 陈洁, 潘雄. 综合半参数核估计与正则化方法的航空重力向下延拓模型分析[J]. 武汉大学学报 ● 信息科学版, 2020, 45(10): 1563-1569. doi: 10.13203/j.whugis20180491
WU Fengfeng, HUANG Haijun, REN Qingyang, FAN Wenyou, CHEN Jie, PAN Xiong. Analysis of Downward Continuation Model of Airborne Gravity Based on Comprehensive Semi-parametric Kernel Estimation and Regularization Method[J]. Geomatics and Information Science of Wuhan University, 2020, 45(10): 1563-1569. doi: 10.13203/j.whugis20180491
Citation: WU Fengfeng, HUANG Haijun, REN Qingyang, FAN Wenyou, CHEN Jie, PAN Xiong. Analysis of Downward Continuation Model of Airborne Gravity Based on Comprehensive Semi-parametric Kernel Estimation and Regularization Method[J]. Geomatics and Information Science of Wuhan University, 2020, 45(10): 1563-1569. doi: 10.13203/j.whugis20180491
  • 航空重力测量与传统的重力测量方法相比,具有测量速度快、范围广、成本低等优势,除此之外,还可在沙漠、沼泽、原始森林等常规重力测量无法作业的地方进行,从而使航空重力成为大规模重力数据采集的主要手段[1-2]

    航空重力测量技术一个重要的环节就是数据的向下延拓,延拓精度的高低将直接影响到重力数据的应用。向下延拓主要采用的是求解逆泊松积分[3]方程,由于求逆不稳定,难以得到准确的延拓解,该问题属于不适定问题。为了减弱病态性的影响,国内外学者提出了很多解决思路,主要有3类:第1类为直接求逆泊松积分方程,以在求逆过程中加入正则化参数抑制噪声为主,如Tikhonov正则化算法(TIK)[4]、截断奇异值分解法[5]等;第2类为间接求逆,如矩谐分析法[6]、最小二乘配置(least squres,LS)[7-9]等,避开了直接对泊松积分求逆,但以LS为例,解算过程中协方差阵仍涉及求逆,存在解不稳定的情况。第3类为非求逆法,如利用卫星测高重力向上延拓和超高阶位模型计算改正数、联合使用超高位模型和地形信息来确定延拓改正数[10]及直接代表法[11]等。除求逆以外,在参数解算过程中,观测数据还受到系统误差的影响,如格值标定误差、加速度改正的误差[12-13]等。为了保证延拓精度,在进行延拓求解时,需要事先分离出系统误差,一般情况是引入高质量的外部数据来改正系统误差。若无高质量的外部数据,上述3类延拓方法都会受系统误差影响,影响延拓的精度。除此之外,在建立逆泊松积分向下延拓模型计算过程中,还存在积分离散化问题[14],若使用不同的离散方法,得出的延拓结果会存在一定的差异,积分离散会对系数矩阵造成干扰,导致存在一定的模型误差[15]

    为了解决上述误差对延拓精度的影响,将重力数据中的系统误差和离散化造成的模型误差用非参数分量表达,在无外部数据的情况下,引入核函数[16],建立基于半参数核估计方法的重力向下延拓模型,为了改善逆泊松积分离散后的设计矩阵的病态影响,引入正则化方法,提出了综合半参数核估计和正则化方法的逆泊松积分延拓方法,在无外部数据时,向下延拓的同时分离系统误差,并与传统的LS法、TIK法等延拓模型进行了精度比较,利用实测数据和模拟数据验证了本文提出的半参数正则化模型向下延拓的可靠性、适用性。

    • 根据泊松积分公式,将地面重力异常数据向上延拓,有:

      $$\text{ }\!\!\Delta\!\!\text{ }{{g}^{\text{*}}}\left( r, \theta , \lambda \right)=\frac{\left( {{R}^{2}}\left( {{r}^{2}}-{{R}^{2}} \right) \right)}{4\text{ }\!\!\pi\!\!\text{ }r}\iint\limits_{\sigma }{\frac{\text{ }\!\!\Delta\!\!\text{ }g\left( R, {\theta }', {\lambda }' \right)}{{{\rho }^{3}}}\text{d}\sigma }$$ (1)

      式(1)将积分面视为球面,近似表达大地水准面,因此需先将地面数据解析延拓至大地水准面上,记为${\rm{\Delta }}g\left( {R, \theta ', \lambda '} \right) $,而${\rm{\Delta }}{g^{\rm{*}}}\left( {r, \theta , \lambda } \right) $则为空中重力异常,$r = R + h $为地球平均半径R与飞行高度h的和;θ'、λ'与θλ分别表示地面点和空中点的纬度、经度;dσ为格网面积,$\rho = \sqrt {{r^2} + {R^2} - 2rR{\rm{cos}}\psi } $为${\rm{\Delta }}{g^{\rm{*}}}\left( {r, \theta , \lambda } \right) $到${\rm{\Delta }}g\left( {R, \theta ', \lambda '} \right) $的距离,Ψ为球面角度,${\rm{cos}}\psi = {\rm{cos}}\theta {\rm{cos}}\theta ' + {\rm{sin}}\theta {\rm{sin}}\theta '{\rm{cos}}\left( {\lambda ' - \lambda } \right) $。

      将泊松积分离散后,矩阵形式的泊松积分可表示为:

      $${\rm{\Delta }}g_i^{\rm{*}} = \mathop \sum \limits_{j = 1}^n {B_{ij}}{\rm{\Delta }}{g_j}$$ (2)

      式中,i=1, 2…m, m为空中重力异常点个数;n为地面重力异常点个数。向下延拓即求式(2)的逆问题。

    • 向下延拓过程中会放大随机噪声、系统误差及逆泊松积分离散化造成的模型误差,为了减弱系统误差、模型误差的影响,可以将其都归到非参数分量,用非参数分量s(ti)表达,此时,式(2)变为:

      $${\rm{\Delta }}g_i^{\rm{*}} = \mathop \sum \limits_{j = 1}^n {B_{ij}}{\rm{\Delta }}{g_j} + s\left( {{t_i}} \right)$$ (3)

      表达为向量形式为:${\mathit{\boldsymbol{g}}^{\rm{*}}} = \mathit{\boldsymbol{BX}} + \mathit{\boldsymbol{s}} $,其中,$ {\mathit{\boldsymbol{g}}^{\rm{*}}} = {[{\rm{\Delta }}g_1^{\rm{*}}, {\rm{\Delta }}g_2^{\rm{*}} \cdots {\rm{\Delta }}g_m^{\rm{*}}]^{\rm{T}}}$为空中重力观测量,$\mathit{\boldsymbol{X}} = {[{\rm{\Delta }}{g_1}, {\rm{\Delta }}{g_2} \cdots {\rm{\Delta }}{g_n}]^{\rm{T}}} $为待求地面重力,$\mathit{\boldsymbol{B}} = {[{b_1}, {b_2} \cdots {b_m}]^{\rm{T}}} $为逆泊松积分离散后的设计矩阵,$\mathit{\boldsymbol{s}}{\rm{ = }}{\left[ {s\left( {{t_1}} \right), s\left( {{t_2}} \right) \cdots s({t_m})} \right]^{\rm{T}}} $是非参数分量。

      求解时,发现当矩阵BTPB病态时,补偿最小二乘法不能得到理想效果,因此为了求解参数的估计值,首先假设X已知,基于$ \{ {t_i}, g_i^{\rm{*}} - \mathit{\boldsymbol{b}}_i^{\rm{T}}\mathit{\boldsymbol{X}}\} _{i = 1}^m$可以对s(tk)做核估计:

      $${{\hat s}_k} = \mathop \sum \limits_{i = 1}^m {W_i}\left( {{t_k}} \right)\left( {g_i^{\rm{*}} - \mathit{\boldsymbol{b}}_i^{\rm{T}}\mathit{\boldsymbol{X}}} \right)$$ (4)

      式中,一个系统误差对应一个${W_i}\left( {{t_k}} \right) $,${W_i}\left( {{t_k}} \right) $为核权函数,是整个样本相对于点tk的权,反映估计${{\hat s}_k} $时观测量作用大小。对选定点tk有:${W_i}\left( {{t_k}} \right) = K\left( {\left( {{t_k} - {t_i}} \right){h^{ - 1}}} \right)/\mathop \sum \limits_{j = 1}^m K\left( {\left( {{t_k} - {t_j}} \right){h^{ - 1}}} \right) $,其中,$ i, j = 1, 2 \cdots m{\rm{}};{\rm{}}k = 1, 2 \cdots m$,m为航空重力数据的个数;K(·)为选取的核函数;h为窗宽参数,且h >0。

      此时,由式(4)可计算残差为:

      $${V_i} = \mathit{\boldsymbol{b}}_i^{\rm{T}}\mathit{\boldsymbol{X}} + {{\hat s}_i} - g_i^{\rm{*}}$$ (5)

      将式(4)代入,并令$W = {({W_i}({t_j}))_{m \times m}} $,则式(5)可以写成如下向量形式:

      $$ \mathit{\boldsymbol{V}} = \left( {\mathit{\boldsymbol{I}} - \mathit{\boldsymbol{W}}} \right)\left( {\mathit{\boldsymbol{BX}} - {\mathit{\boldsymbol{g}}^{\rm{*}}}} \right) $$ (6)

      令$\mathit{\boldsymbol{\tilde B}} = \left( {\mathit{\boldsymbol{I}} - \mathit{\boldsymbol{W}}} \right)\mathit{\boldsymbol{B}} $,${{\mathit{\boldsymbol{\tilde g}}}^{\rm{*}}} = \left( {\mathit{\boldsymbol{I}} - \mathit{\boldsymbol{W}}} \right){\mathit{\boldsymbol{g}}^{\rm{*}}} $,由于${({{\mathit{\boldsymbol{\tilde B}}}^{\rm{T}}}\mathit{\boldsymbol{P\tilde B}})^{ - 1}} $可能会存在奇异,为此,引入正则化方法减弱病态,具体为:

      $$ {(\mathit{\boldsymbol{\tilde B\hat X}} - {{\mathit{\boldsymbol{\tilde g}}}^{\rm{*}}})^{\rm{T}}}\mathit{\boldsymbol{P}}\left( {\mathit{\boldsymbol{\tilde B\hat X}} - {{\mathit{\boldsymbol{\tilde g}}}^{\rm{*}}}} \right) + \beta {{\mathit{\boldsymbol{\hat X}}}^{\rm{T}}}\mathit{\boldsymbol{\hat X}} = {\rm{min}} $$ (7)

      式中,β为正则化参数。构造拉格朗日函数,对参数分别求解,得到参数估计值:

      $$ \mathit{\boldsymbol{\hat X}} = {({{\mathit{\boldsymbol{\tilde B}}}^{\rm{T}}}\mathit{\boldsymbol{P\tilde B}} + \beta \mathit{\boldsymbol{I}})^{ - 1}}{{\mathit{\boldsymbol{\tilde B}}}^{\rm{T}}}P{{\mathit{\boldsymbol{\tilde g}}}^{\rm{*}}} $$ (8)

      将式(8)代入式(4)中,即可得到非参数分量${\hat s} $的估计值:

      $$ \mathit{\boldsymbol{\hat s}} = \mathit{\boldsymbol{W}}\left( {{\mathit{\boldsymbol{g}}^{\rm{*}}} - \mathit{\boldsymbol{B\hat X}}} \right) $$ (9)

      从式(8)、式(9)可以看出,估计值精度的高低与窗宽参数和正则化参数的合理选择有关,参数的选择使用广义交叉核实法(generalized cross-validation,GCV)来选取。GCV函数定义为:

      $$ {\rm{GCV}}\left( {h, \beta } \right) = \frac{{{\mathit{\boldsymbol{V}}^{\rm{T}}}\mathit{\boldsymbol{PV}}}}{{{{(1 - {\rm{tr}}(\mathit{\boldsymbol{H}}\left( {h, \beta } \right))/m)}^2}}} $$ (10)

      式中,H(h, β)为影响矩阵,由$\mathit{\boldsymbol{BX}} + \mathit{\boldsymbol{s}} = \mathit{\boldsymbol{H}}\left( {h, \beta } \right){\mathit{\boldsymbol{g}}^{\rm{*}}} $定义;m为观测值个数;tr(H)为矩阵H的迹。H(h, β)的定义如下:

      $$ \begin{array}{*{20}{l}} {\mathit{\boldsymbol{H}}\left( {h, \beta } \right) = \mathit{\boldsymbol{\tilde B}}{{({{\mathit{\boldsymbol{\tilde B}}}^{\rm{T}}}\mathit{\boldsymbol{P\tilde B}} + \mathit{\boldsymbol{\beta I}})}^{ - 1}}{{\mathit{\boldsymbol{\tilde B}}}^{\rm{T}}}(\mathit{\boldsymbol{I}} - {\rm{}}}\\ {{\rm{}}\mathit{\boldsymbol{W}}\left( h \right)) + \mathit{\boldsymbol{W}}\left( h \right)} \end{array} $$ (11)

      两个参数的求解方法是设置GCV双循环,求出GCV函数值最小时对应的窗宽和正则化参数。

    • 为了验证本文提出的半参数核估计正则化向下延拓方法的延拓精度与可靠性,设计了一个仿真算例和一个实测算例,分别以LS法、TIK法和本文方法进行向下延拓,并进行对比分析。

    • 数据选择美国南部一个4°×1°区域作为实验区,其区域为27°N~31°N、87°W~88°W,使用EGM2008(earth gravity model 2008)模型2 160阶计算了该区域飞行高度为3.3 km与6.3 km处分辨率为3′×3′的空中重力异常,为评估延拓效果,同时计算了零高度面的重力异常作为基准值。表 1统计了仿真区域EGM2008模拟的重力异常特性。

      表 1  仿真区域各项重力异常特性的统计

      Table 1.  Characteristics of Gravity Anomalies at Simulation Area

      高度/km 最小值
      /mGal
      最大值
      /mGal
      平均值
      /mGal
      标准差
      /mGal
      0 -104.35 35.67 -29.78 32.93
      3.3 -90.93 28.03 -28.93 30.07
      6.3 -82.97 22.73 -28.19 27.84
    • 1)首先向空中重力加入系统误差及随机噪声生成仿真空中重力异常值,利用移去-恢复方法,取EGM2008重力场模型的360阶作为参考模型,得到参考重力异常值,从仿真空中重力异常值中移去,得到残差空中重力异常值。

      2)分别采用LS法、TIK法与半参数估计模型(semiparametric estimation model,SEM)将残差空中重力异常延拓到零高度面并恢复参考场,将其同基准值进行比较,从而获取相应的计算模型精度评估参数。

      3)为评估半参数模型在无外部数据时分离系统误差的能力,比较模拟系统误差与半参数模型分离的系统误差的大小及精度。

    • 为了探究随机噪声及不同类型系统误差对延拓结果的影响,设计了以下两种方案生成系统误差。

      方案1   生成以线性项为主要变化量的系统误差Simu1,其常数项为1 mGal,漂移系数为0.002 5 mGal/s,周期项(正弦函数)振幅为1.5 mGal、频率为0.001 Hz;

      方案2   生成以周期项为主要变化量的系统误差Simu2,其常数项为3 mGal,漂移系数为-0.001 mGal/s,周期项(正弦函数)振幅为5 mGal,频率为0.003 Hz。

      向空中重力异常加入系统误差及随机噪声仿真,仿真数据分别用LS、TIK、SEM按计算步骤向下延拓,3种方法延拓值与基准值差值的具体统计信息见表 2

      表 2  不同方法延拓值与基准值差值

      Table 2.  Differences Between the Results and Real Data from Different Models

      延拓高度
      /km
      系统误差 随机噪声
      /km
      LS/mGal TIK/mGal SEM/mGal
      平均值 均方根 平均值 均方根 平均值 均方根
      3.3 Simu1 ±1 3.38 5.01 3.06 3.96 0.06 2.23
      ±3 3.49 10.77 2.52 4.66 0.08 3.92
      6.3 Simu1 ±1 3.49 9.65 3.20 4.41 0.06 2.78
      ±3 3.60 26.28 2.83 5.15 0.08 4.26
      3.3 Simu2 ±1 2.42 5.75 1.94 4.46 0.10 3.25
      ±3 2.47 11.02 1.70 4.95 0.10 4.35
      6.3 Simu2 ±1 2.50 9.97 2.13 5.15 0.11 3.76
      ±3 2.55 25.98 1.87 5.53 0.11 4.79
      注:Simu1、Simu2分别为方案1、方案2中的模拟系统误差

      向下延拓属不适定问题,法矩阵一般有较强的病态特征,从表 2可看出,LS法由于病态原因,随机噪声被放大,其精度降低。TIK和SEM表现出了较强的抗干扰能力。在延拓高度相同情况下,使随机噪声由±1 mGal增加至±3 mGal,系统误差为Simu1时,TIK和SEM延拓精度平均降低0.7 mGal和1.6 mGal左右;系统误差为Simu2时,TIK和SEM延拓精度平均降低0.4 mGal和1 mGal左右。在误差相同情况下,延拓高度从3.3 km增加至6.3 km,TIK精度降低0.5~0.6 mGal左右,SEM精度降低0.3~0.4 mGal,可知TIK和SEM能改善病态问题,对随机噪声都有一定的抑制作用。

      从受系统误差影响大小来看,只有SEM差值平均值一直在零附近,明显优于LS与TIK,当延拓受Simu1影响时,LS和TIK差值平均值分别偏移至3.48 mGal和2.8 mGal,当延拓受Simu2影响时,两方法差值平均值分别偏移2.48 mGal和1.9 mGal,可知LS和TIK受到了系统误差影响,降低了延拓精度。当系统误差分别为Simu1与Simu2时,SEM平均延拓精度分别高于正则化方法1.25 mGal、0.93 mGal。从差值的平均值与延拓精度两方面来看,说明TIK配合移去-恢复法对系统误差稍有改善,但SEM更加有效地分离了系统误差。

      为了使误差分布情况更直观,以表 2中两种方案在高度6.3 km及随机噪声±3 mGal情况下为例绘制误差分布图,由于LS延拓结果在此条件下已失真,以下只绘制TIK与SEM模型的误差分布图,如图 1所示。

      图  1  向下延拓值与基准值的差值统计

      Figure 1.  Differences Between the Results of Downward Continuation and Real Data

      对比图 1方案1,TIK虽改善了病态,延拓结果较好,但能发现延拓结果受到了系统误差中线性项的影响,随纬度从31°N~27°N,其延拓差值愈发增大;对比图 1方案2,当系统误差变化以周期项为主时,TIK的延拓差值大小在纬度上呈带状分布,两种差值分布情况与加入的系统误差相符合。从图 1中可发现系统误差对于半参数模型的延拓结果来说影响并不明显。表 3进一步列出了利用SEM估计观测数据中的系统误差情况。

      表 3  利用SEM估计观测数据中的系统误差情况

      Table 3.  Statistics of Systematic Errors Estimated by SEM Data

      延拓高度/km 系统
      误差
      随机噪声
      /mGal
      最小值
      /mGal
      最大值
      /mGal
      平均值
      /mGal
      标准差
      /mGal
      Simu1 - 1.01 5.68 3.27 1.47
      3.3 ±1 1.15 5.48 3.22 1.35
      ±3 1.37 5.51 3.31 1.32
      6.3 ±1 1.23 5.29 3.22 1.22
      ±3 1.30 5.47 3.32 1.27
      Simu2 - -3.38 7.88 2.37 3.61
      3.3 ±1 -3.57 8.04 2.32 3.14
      ±3 -3.53 7.82 2.49 3.24
      6.3 ±1 -3.38 7.98 2.32 3.14
      ±3 -3.32 7.75 2.50 3.23
      注:Simu1、Simu2分别为方案1、方案2中的模拟系统误差

      为直观表现分离效果,同样以上述两种方案在高度6.3 km及随机噪声±3 mGal情况下为例绘制系统误差估计如图 2所示。

      表 3中可以看出,方案1中,在3.3 km与6.3 km延拓高度下,系统误差估计值与真值的标准差相差不过0.2 mGal;而在方案2中,其标准差相差在0.5 mGal左右,SEM分离的系统误差能在一定程度上吻合真值。但可以发现,系统误差分离效果随随机噪声的增大而降低。图 2反映出SEM在面对不同种类系统误差时,都有较为可靠的系统误差分离能力。

      图  2  两方案6.3 km向下延拓系统误差估计情况

      Figure 2.  Systematic Errors Estimation of 6.3 km Downward Continuation of Two Schemes

    • 考虑到仿真数据算例皆在无外部数据情况进行,以实测数据来进一步验证SEM分离系统误差、提高延拓精度的有效性。首先不引入外部数据,将实测数据以本文方法进行向下延拓,并与LS及TIK比较分析;然后再引入外部数据,并用不同方法进行系统误差改正;最后将改正后的数据用TIK向下延拓进行对比分析。

      ES10项目数据由美国国家大地测量局于2016年6月发布,位于美国卡罗来纳州的山区,有较大的地形起伏。其发布数据主测线45条,交叉测线5条,区域为34°N~37°N、80°W~84°W,飞行高度平均在5.4 km。本文目前缺少该地区实测地面重力异常,采用重力模型生成的零高度面重力异常代替地面实测值。

    • 1)在空中实测数据未接受系统误差改正的情况下,将数据格网化为分辨率5′×5′的网格空中重力异常数据,使用EGM2008模型2 160阶计算同范围、同分辨率的零高度重力异常作为参考值。分别采用LS、TIK与SEM将空中重力异常向下延拓后与基准值进行比较。

      2)利用SEM获得延拓解的同时分离出了系统误差,以此系统误差改正空中重力异常。

      3)分别选用交叉点平差两步法、整条测线改正和有外部数据时的SEM来进行系统误差估计。为对比系统误差改正效果,将步骤2)中改正后的数据与经过3种方法改正后的空中实测数据进行向下延拓实验。

    • 首先分析在无外部数据情况下的延拓结果,从表 4可看到,就差值平均值而言,SEM受系统误差影响最小,TIK较LS稍有改善。就差值均方根而言,半参数模型比另两种模型分别高出0.64 mGal、0.33 mGal,这与模拟算例结果是相符的。

      表 4  无外部数据延拓值与基准值差值/mGal

      Table 4.  Differences Between the Continuation Results and Real Data from Different Models (No External Data)/mGal

      模型 最小值 最大值 平均值 均方根
      LS -29.39 33.40 2.01 6.74
      TIK -27.69 31.54 -1.42 6.43
      SEM -29.67 28.81 -0.05 6.10

      为进一步说明SEM无外部数据也能有效分离系统误差,引入外部数据进行改正,将改正后的数据向下延拓,并与本文方法对比。首先采用交叉点两步平差法来进行测线系统误差改正,表 5给出了平差前后交叉点不符值的相关信息统计。

      表 5  平差前后交叉点不符值统计/mGal

      Table 5.  Dis crepancy Statistics of Cross Points Before and After Adjustment/mGal

      平差前后 最小值 最大值 平均值 标准差 均方根
      平差前 -6.43 5.56 -0.33 2.09 2.12
      平差后 -4.46 5.11 -0.002 1.47 1.47

      然后采用整条测线改正,EGM2008模型建模时加入了大量北美地区重力数据,因此EGM2008在北美地区有较高精度,便直接用高阶次重力场来进行测线系统误差改正。表 6给出了部分主测线与EGM2008模型的偏差。

      表 6  部分主测线与EGM2008模型偏差统计/mGal

      Table 6.  Deviation Statistics of Part Main Survey Lines and EGM2008/mGal

      测线 偏差 测线 偏差 测线 偏差
      101 3.11 108 1.54 113 2.39
      102 2.19 109 2.46 114 1.21
      104 -0.02 110 2.46 115 2.10
      105 2.81 111 3.13 116 2.84
      106 1.76 112 1.17 117 2.53

      半参数估计模型在无外部数据时,用核估计正则化方法分离系统误差,在有外部数据时,其能更准确地分离系统误差,将实测空中重力数据减去外部数据,以此代替实测值,向下延拓来分离系统误差。分别将无外部数据时半参数估计模型改正、交叉点平差改正、整条测线改正与有外部数据时半参数估计模型改正后所对应的数据向下延拓,此处为方便对比延拓精度,均只采用Tikhonov正则化方法来进行延拓,延拓结果见表 7

      表 7  系统误差改正后延拓值与基准值差值/mGal

      Table 7.  Difference Between the Continuation Results and the Reference Values After the Systematic Error Correction/mGal

      延拓数据 最小值 最大值 平均值 均方根
      无外部数据
      SEM改正
      -29.71 30.34 -0.05 6.10
      交叉点平
      差改正
      -27.23 30.34 1.79 6.37
      整条测试改正 -29.92 29.31 -0.26 6.12
      有外部数据
      SEM改正
      -29.34 28.38 -0.05 6.08

      未改正前,正则化方法延拓精度为6.43 mGal,对比表 7可知,交叉点平差后的数据延拓精度最低,且延拓结果中仍存在系统误差,查阅资料发现ES10交叉测线精度低于正常测线,交叉点平差后导致正常测线数据受到了污染。

      数据改正后,正则化方法延拓精度均提升0.1~0.3 mGal。就差值平均值而言,SEM改正不论是否存在外部数据,其受系统误差的影响都较小,而整条测线改正由于使用了高精度的外部数据,其受系统误差的影响也较小;就差值均方根而言,无外部数据的SEM改正效果与整条测线改正的效果相当,有外部数据的SEM改正后搭配正则化方法延拓精度最高,为6.08 mGal。

      但要补充的是,未经处理的数据使用SEM改正延拓精度也达到了6.10 mGal,无外部数据的SEM向下延拓结果与引入外部数据改正后的正则化方法向下延拓结果精度相当。实测数据向下延拓实验表明,本文方法在无外部数据情况下,能够较为有效地分离系统误差达到更高的精度。

    • 通过仿真数据与实测空中重力数据计算结果对比分析,有以下结论:

      1)高度越高,延拓精度越低,线性项为主的系统误差的延拓值优于周期项为主的系统误差的延拓值。

      2)仿真模拟实验表明,SEM与TIK都可减弱随机噪声的影响,对比延拓差值的平均值和均方根,SEM能更有效地分离出系统误差,提高延拓精度。

      3)通过实测空中重力数据实验,在无外部数据的情况下,SEM无需事先处理即可直接向下延拓,其延拓精度与引入外部数据事先改正系统误差后的TIK延拓精度相当;而在有外部数据时,SEM分离系统误差更准确,结合TIK,延拓精度达到6.08 mGal。

参考文献 (16)

目录

    /

    返回文章
    返回