CN113378471B - 一种顾及重力地形改正的gps高程拟合方法及*** - Google Patents

一种顾及重力地形改正的gps高程拟合方法及*** Download PDF

Info

Publication number
CN113378471B
CN113378471B CN202110692845.1A CN202110692845A CN113378471B CN 113378471 B CN113378471 B CN 113378471B CN 202110692845 A CN202110692845 A CN 202110692845A CN 113378471 B CN113378471 B CN 113378471B
Authority
CN
China
Prior art keywords
terrain
elevation
correction
gps
gravity
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202110692845.1A
Other languages
English (en)
Other versions
CN113378471A (zh
Inventor
钟波
谭江涛
李贤炮
刘滔
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Wuhan University WHU
Original Assignee
Wuhan University WHU
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Wuhan University WHU filed Critical Wuhan University WHU
Priority to CN202110692845.1A priority Critical patent/CN113378471B/zh
Publication of CN113378471A publication Critical patent/CN113378471A/zh
Application granted granted Critical
Publication of CN113378471B publication Critical patent/CN113378471B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G06F30/27—Design optimisation, verification or simulation using machine learning, e.g. artificial intelligence, neural networks, support vector machines [SVM] or training a model
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01C—MEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
    • G01C5/00—Measuring height; Measuring distances transverse to line of sight; Levelling between separated points; Surveyors' levels
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
    • G06N3/00—Computing arrangements based on biological models
    • G06N3/02—Neural networks
    • G06N3/04—Architecture, e.g. interconnection topology
    • G06N3/045—Combinations of networks
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
    • G06N3/00—Computing arrangements based on biological models
    • G06N3/02—Neural networks
    • G06N3/08—Learning methods
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00—Details relating to CAD techniques
    • G06F2111/10—Numerical modelling
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
    • G06F2119/14—Force analysis or force optimisation, e.g. static or dynamic forces

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Evolutionary Computation (AREA)
  • General Physics & Mathematics (AREA)
  • Software Systems (AREA)
  • Artificial Intelligence (AREA)
  • General Engineering & Computer Science (AREA)
  • Computing Systems (AREA)
  • Health & Medical Sciences (AREA)
  • General Health & Medical Sciences (AREA)
  • Molecular Biology (AREA)
  • Computational Linguistics (AREA)
  • Biophysics (AREA)
  • Biomedical Technology (AREA)
  • Mathematical Physics (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Data Mining & Analysis (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Medical Informatics (AREA)
  • Computer Hardware Design (AREA)
  • Geometry (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Position Fixing By Use Of Radio Waves (AREA)

Abstract

本发明提供一种顾及重力地形改正的GPS高程拟合方法，用于山区地形，包括根据测区的数字高程模型计算格网点处布格重力异常的地形改正和层间改正；将布格重力异常的地形改正和层间改正相加，转换得到格网点处高程异常中地形因素的影响；进行重力地形改正扣除GPS水准点处高程异常中地形因素的影响；对地形改正后GPS水准点的高程异常进行数值拟合，得到地形改正后的局部似大地水准面模型；将地形改正后的局部似大地水准面模型与对应位置计算得到的高程异常中的地形影响结合，得到真实的局部似大地水准面，作为最终所得顾及重力地形改正的GPS高程拟合成果。本发明数据获取成本较重力测量方法低，可以有效提高山区的GPS高程拟合精度。

Description

一种顾及重力地形改正的GPS高程拟合方法及***
技术领域
本发明涉及一种高精度似大地水准面模型构建方法，特别涉及一种在山地区域考虑地形起伏影响进行重力地形改正的GPS高程拟合方案，属于大地测量与测量工程领域。
背景技术
GPS测量技术所得测点的高程为基于参考椭球面的大地高，而我国采用高程***为基于似大地水准面的正常高***，两者存在一定的差距，若要将GPS测量的高程结果应用于实际工程中，需得到测区的似大地水准面模型，从而把大地高转换为正常高。GPS高程拟合方法是构建局部似大地水准面模型方法中的一种，其主要原理为利用区域中一些已知大地高和正常高的点位数据(称这些点为GPS水准点)，通过数值拟合的方式拟合构建该区域的似大地水准面模型。目前GPS高程拟合的方法有很多种，如多面函数法[1]、Shepard曲面法[2]、多项式拟合法[3]和BP神经网络方法[4]。由于局部区域的高程异常受地形起伏的影响较大，实际应用中常规GPS高程拟合方法在山地区域很难得到高精度的拟合模型，因此，若想获得高精度的区域GPS高程拟合结果，必须考虑地形起伏对高程异常的影响。
重力地形改正[5]主要用于局部重力场的正常重力值求解中，而局部重力场与高程异常之间存在一定的数学关系，依据该数学关系可以构建高精度的局部似大地水准面模型，同时，重力场理论中重力地形改正可以解决地形因素在重力场中的影响。目前，重力地形改正通常采用Stokes方法[6]，其基本原理为，“移去”陆地地形，再将移去后的地形质量以合理的方式“补偿”到似大地水准面上，在地球总质量和质心不变的基础上，实现对地形质量进行调整。调整过程中，移去或补偿的地形质量可以等效为地形引起的重力变化，即重力地形改正项，通过测定区域内测点处的重力值并结合重力改正项即可得到正常重力值，进而得到区域似大地水准面模型。虽然Stokes方法构建的区域似大地水准面模型精度高且可靠性好，但该方法需要精确测量区域地面重力值，数据采集成本较高且耗时耗力。另一种重力地形改正方法则是从布格重力异常的角度进行改正[7]，其基本原理为，将给定范围内的高于测点的地形内的物质“削平”，同时将低于测点的地形用同一密度物质“填满”，从而使得测点处的重力不受周围地形起伏的影响，该步骤称为布格重力异常的地形改正，之后再根据测点所处高程将所有测点归算至统一的基准面上，即改正测点正下方的地形影响，该步骤称为布格重力异常的层间改正。利用布格重力异常来确定区域似大地水准面虽然也需进行重力测量，但相比于Stokes方法，该方法可以直接通过区域高程分布情况以及地质密度改正绝大部分地形起伏的影响。
在山地区域的似大地水准面模型构建中，通常使用“移去-恢复”的组合拟合方法，先对高程异常变化中变化频率较小的部分采用合适的方法拟合，再用另一种方法拟合变化频率较大的残差项，如基于二次曲面的Shepard拟合法[8]或者基于二次曲面的多面函数法[9]，也有使用重力测量与数值拟合相结合的方法，但前者对地形起伏仅进行了拟合而不是改正，精度稍弱，而后者精度虽高可成本也高。因此，该领域亟待提出精度更佳且无需重力测量的似大地水准面模型构建方法。
本发明中涉及的参考文献如下：
[1]Hardy R L.Multiquadric equations of topography and other irregularsurfaces[J].Journal of Geophysical Research,1971,76(8):1905-1915.
[2]丁海勇，杨延珍，王新志，等.基于移动曲面拟合方法的GPS高程转换[J].大地测量与地球动力学，2010，30(6)：82-85.
[3]杨德明.GPS水准多项式曲面拟合模型研究[J].测绘与空间地理信息，2006，29(3)：34-36.
[4]杨明清，靳蕃，朱达成，等.用神经网络方法转换GPS高程[J].测绘学报，1999，28(4)：301-307.
[5]李振海，李琼，林旭.重力地形改正的计算模型研究[J].测绘工程，2011，20(2)：24-26.
[6]STOKES G G.On the variation of gravity at the surface of the Earth[J].Transcation of the Cambridge Philosophical Society,1849(8):672-695.
[7]高庆余.布格重力异常地形校正方法在微重力勘探中的缺点及其纠正方法[J].物探化探计算技术，1998，20(4)：319-327.
[8]侯本军，钟波，王文庆.基于二次曲面的Shepard拟合法在GPS水准中的应用[J].测绘工程，2007，16(1)：36-38.
[9]朱开银，秦岩宾，何友福.二次-多面函数模型的GPS高程拟合及精度分析[J].测绘与空间地理信息，2019，42(2)：117-119.
发明内容
针对重力测量方法成本较高而常规GPS拟合方法山区精度较差的问题，本发明提出了一种顾及重力地形改正的GPS高程拟合方案，使用成本较低的纯数值拟合方法建立局部似大地水准面模型，适用于山区地形。
为了实现上述目的，本发明提出的技术方案为一种顾及重力地形改正的GPS高程拟合方法，用于山区地形的GPS高程拟合，包括步骤如下：
步骤1，根据测区的数字高程模型中包含的大地高信息以及测区格网信息计算格网点处布格重力异常的地形改正TC项和层间改正BP项；
步骤2，将布格重力异常的地形改正和层间改正相加，所得结果δg转换为对似大地水准面的影响，得到格网点处高程异常中地形因素的影响ΔNP；
步骤3，选取测区中的GPS水准点，利用测点上的大地高H和正常高h相减高程异常ξ，再去除GPS水准点处地形因素的影响ΔN，得到经过地形改正后高程异常ξ′；
步骤4，对地形改正后GPS水准点的高程异常进行数值拟合，得到地形改正后的局部似大地水准面模型N'；
步骤5，将步骤4所得地形改正后的局部似大地水准面模型N'与对应位置计算得到的高程异常中的地形影响ΔNP结合，得到真实的局部似大地水准面N，作为最终所得顾及重力地形改正的GPS高程拟合成果。
而且，步骤1中，根据地形起伏与地面点间的关系，直接给出地形改正TC项计算公式如下，
Figure BDA0003127370690000031
其中，G为地心引力常数，ρ(x,y,z)为流动点(x,y,z)的地壳密度，h为流动点的高程，E代表积分区域，r()为距离函数；
层间改正BP项的计算公式为BP(i,j)＝-2πGρ0hij，其中，ρ0为地壳密度。
而且，步骤2中，根据重力异常与高程异常间的关系，某地面点P上布格重力异常的层间改正和地形改正对似大地水准面的影响用Stokes积分公式表示成如下形式，
Figure BDA0003127370690000032
其中，γ为平均正常重力，E代表积分区域，R为地球平均半径，
Figure BDA0003127370690000041
为Stokes核函数，
Figure BDA0003127370690000042
和λP为地面点P的纬度和经度，
Figure BDA0003127370690000043
和λ为点P周围流动点的纬度和经度。
或者，步骤2中，布格重力异常的地形改正和层间改正对高程异常间的关系表示为球谐展开式的形式，
Figure BDA0003127370690000044
其中，G为地心引力常数，γ为平均正常重力，R为地球平均半径，ρ为地壳密度，
Figure BDA0003127370690000045
h为流动点
Figure BDA0003127370690000046
的高程；Ynm为n阶m次的完全规格化的球谐函数；
Figure BDA0003127370690000047
和λP为地面点P的纬度和经度，
Figure BDA0003127370690000048
和λ为点P周围流动点的纬度和经度。
而且，步骤3中，根据格网点处高程异常中地形因素的影响
Figure BDA0003127370690000049
以反距离加权插值的方式得到GPS水准点处高程异常中地形因素的影响ΔN，反距离加权插值表示成如下形式，
Figure BDA00031273706900000410
其中，Pi,i＝1,2,3,4为距离该GPS水准点最近的4个格网点，di为GPS水准点至格网点的平面距离。
而且，步骤4中，对地形改正后的GPS水准点高程异常进行数值拟合采用多项式方法、Shepard方法、多面函数方法或BP神经网络方法。
而且，将步骤4经过数值拟合构建得到的地形改正后的局部似大地水准面模型N'用格网数据形式表达，使得拟合模型的格网点与步骤1所述测区格网一一对应，则在步骤5恢复真实的局部似大地水准面模型N时，直接利用步骤2中的格网点处高程异常中地形因素的影响ΔNP，无需再计算对应位置布格重力异常中的层间改正与地形改正并用Stokes积分转换为高程异常中地形因素的影响。
另一方面，本发明提供一种顾及重力地形改正的GPS高程拟合***，用于实现如上所述的一种顾及重力地形改正的GPS高程拟合方法。
而且，包括以下模块，
第一模块，用于根据测区的数字高程模型中包含的大地高信息以及测区格网信息计算格网点处布格重力异常的地形改正TC项和层间改正BP项；
第二模块，用于将布格重力异常的地形改正和层间改正相加，所得结果δg使用Stokes积分公式转换为对似大地水准面的影响，得到格网点处高程异常中地形因素的影响ΔNP；
第三模块，用于选取测区中的GPS水准点，利用测点上的大地高H和正常高h相减高程异常ξ，再去除GPS水准点处地形因素的影响ΔN，得到经过地形改正后高程异常ξ′；
第四模块，用于对地形改正后GPS水准点的高程异常进行数值拟合，得到地形改正后的局部似大地水准面模型N'；
第五模块，用于将第四模块所得地形改正后的局部似大地水准面模型N'与对应位置计算得到的高程异常中的地形影响ΔNP结合，得到真实的局部似大地水准面N，作为最终所得顾及重力地形改正的GPS高程拟合成果。
或者，包括处理器和存储器，存储器用于存储程序指令，处理器用于调用存储器中的存储指令执行如上所述的一种顾及重力地形改正的GPS高程拟合方法。
与现有技术相比本发明具有以下优点和有益效果：
1、基于纯数值拟合，数据获取成本较重力测量方法低。
当前局部似大地水准面模型的构建方法有重力测量方法和GPS高程拟合方法，重力测量方法精度虽高但测量成本也高，而GPS高程拟合方法测量成本较低但应用于山区时精度较差。本发明将重力测量理论中关于重力地形改正方面的计算方法，结合常规GPS高程拟合方法，只需成本较低的GPS测量和水准测量来获取数据，节省了测量所需的人力和仪器成本。
2、方案可以有效提高山区的GPS高程拟合精度。
常规的GPS高程拟合方法有对不同模型在山区适用情况的研究，但并未从地形因素影响的原理上予以解决，得到的改进效果不尽人意。本发明从重力值与高程异常的关系出发，将重力地形改正转换高程异常中地形因素的改正，进而用于山区的GPS高程拟合，实际应用情况也表明，对于不同GPS拟合方法均可以有效提高在山区的GPS高程拟合精度(改进效果见图4)。
将本发明方法应用于实例测区并与未经改正的GPS拟合方法的检核点拟合残差进行比较，本发明方法所得拟合模型的外符合精度有约36％～52％的改进效果，验证了本发明方法的有效性。
本发明方案实施简单方便，实用性强，解决了相关技术存在的实用性低及实际应用不便的问题，能够提高用户体验，具有重要的市场价值。
附图说明
图1为本发明实施例的方法流程图；
图2为本发明实施例的GPS水准点分布情况以及使用的数字高程模型；
图3为本发明实施例布格重力异常的地形改正和层间改正原理示意图；
图4为本发明实施例提供的顾及重力地形改正的方法和未经改正的常规GPS高程拟合方法各检核点的残差比较情况。
具体实施方式
以下结合附图和实施例具体说明本发明的技术方案。
本发明提供了一种顾及重力地形改正的GPS高程拟合方法，其特点为将布格重力异常中的地形改正和层间改正项通过Stokes积分公式转换为对似大地水准面的影响，实现对GPS水准数据中地形起伏的改正，使其改正后的似大地水准面更平滑从而使GPS高程拟合方法的拟合效果更好。
图1给出了本发明方法的流程，结合图2所示测区的数据进行实验，验证本发明的优化效果。实验测区面积约11500km2，地形起伏落差4000m以上，测区中部有一条山脉，与周围地形的海拔差接近1000m，共有144个GPS水准点，设置其中的22个为检核点来检验所得的区域似大地水准面模型结果。
参见图1，本发明实施例提供一种顾及重力地形改正的GPS高程拟合方法，使用成本较低的纯数值拟合方法建立局部似大地水准面模型，适用于山区地形，实现过程包括以下步骤：
步骤1，根据测区的数字高程模型中包含的大地高信息以及设计的测区格网信息计算格网点处布格重力异常的地形改正TC项和层间改正BP项。
实施例根据测区点位分布，设计测区格网的分辨率为3′×3′，结合图2所示数字高程模型计算布格重力异常的地形改正和层间改正。
如图3所示，方格网填充的区域为地形改正的部分，密集点填充的区域为层间改正的部分。地形改正的基本原理为，将研究点P点周围区域(地面)高于P点所在平面(布格面)的A区域地形质量ΔmTC-移除，将低于P点所在平面亏损的B区域地形质量ΔmTC+补上，因填充或移除地形质量对P点重力值的影响即为地形改正TCP，根据地形起伏与地面点间的关系，直接给出P点地形改正项计算公式如下，
Figure BDA0003127370690000071
其中，G为地心引力常数，取值为3.98603×1014m3/s2，(x,y,z)为A区域和B区域内流动点空间坐标，ρ(x,y,z)为流动点的地壳密度，实施例采用通常取值2.67g/cm3，h为流动点处地形的高程。(xP,yP)为P点平面坐标，hP为P点及其在大地水准面上的投影点P0之间的距离，即P点高程，单位为m，E代表积分区域，即A区域和B区域平面范围之和，
r()为距离函数，
Figure BDA0003127370690000072
最终求得的TCP单位为mGal。
在层间改正中，称P点所在平面为布格面，经过地形改正后其与大地水准面平行，两者之间的地形质量称为布格片，将P点周围区域布格面与大地水准面之间的地形质量ΔmBP移除，因布格片移除对P点重力值的影响即为层间改正，其具体的计算方法为BPP＝-2πGρhP，其中，ρ为地壳密度，实施例采用通常取值2.67g/cm3，由于层间改正计算公式中只有hP为变量，可以将其简化为BPP＝-0.1119hP，hP的单位为m。最终求得的BPP单位为mGal。
步骤2，将布格重力异常的地形改正和层间改正相加得到结果δg＝BP+TC，相加结果使用Stokes积分公式转换为对似大地水准面的影响，即格网点处高程异常中地形因素的影响ΔNP。
实施例中，根据重力异常与高程异常间的关系，某地面点P上布格重力异常的层间改正和地形改正对似大地水准面的影响用Stokes积分公式可以表示成具体转换公式为，
Figure BDA0003127370690000073
其中，
Figure BDA0003127370690000074
和λP为P点的纬度和经度，
Figure BDA0003127370690000075
和λ为点P周围流动点的纬度和经度，E代表积分区域，γ为平均正常重力，取值为9.80m/s2，R为地球平均半径，取值为6371km，
Figure BDA0003127370690000076
为Stokes核函数，具体实施时可应用严密的球面核函数和一维FFT技术计算。
步骤2中，布格重力异常的地形改正和层间改正对高程异常间的关系还可以表示为球谐展开式的形式，
Figure BDA0003127370690000081
其中，G为地心引力常数，ρ为地壳密度，实施例采用通常取值2.67g/cm3；
Figure BDA0003127370690000082
h为流动点
Figure BDA0003127370690000083
的高程；Ynm为n阶m次的完全规格化的球谐函数。
步骤3，选取测区中的GPS水准点，利用测点上的大地高H和正常高h相减高程异常ξ，再扣除GPS水准点处地形因素的影响ΔN，得到经过地形改正后高程异常ξ′。
选取测区中的GPS水准点，利用测点上的大地高H和正常高h获得高程异常ξ，ξ＝H-h。格网点Pi处高程异常中地形因素的影响
Figure BDA0003127370690000084
以反距离加权插值的方式可以得到GPS水准点处高程异常中地形因素的影响ΔN，反距离加权插值可以表示成如下形式，
Figure BDA0003127370690000085
其中，Pi(i＝1,2,3,4)为距离该GPS水准点最近的4个格网点，di为GPS水准点至格网点的平面距离。利用求得的GPS水准点处高程异常中地形因素的影响对所有GPS水准点的高程异常做地形改正，ξ′＝ξ-ΔN，得到经过地形改正后高程异常ξ′。
步骤4，对地形改正后的GPS水准点高程异常进行数值拟合，构建经地形改正后的局部似大地水准面模型N'。
步骤4中，由于已经考虑了重力地形改正，后续的数值拟合方法可以采用多项式方法、Shepard方法、多面函数方法、BP神经网络方法等适用于平缓曲面拟合的方法。
本实施例选取了三种GPS高程拟合方法，分别为多面函数法、基于二次曲面的Shepard方法和基于二次曲面的多面函数方法，从而检验重力地形改正对于不同GPS高程拟合方法的适用性，证明本发明方法的可靠性。
步骤5，将步骤4得到地形改正后的局部似大地水准面模型N'与步骤2中格网点处高程异常中地形因素的影响ΔNP相加，得到真实的局部似大地水准面模型N，N＝N′+ΔNP，即为本发明方法最终得到的顾及重力地形改正的GPS高程拟合结果。
步骤4中，经过数值拟合构建得到的地形改正后的局部似大地水准面模型N'可以用格网数据形式表达，使得拟合模型的格网点与步骤1设计的格网一一对应，则在步骤5恢复真实的局部似大地水准面模型N时，可以直接利用步骤2中的格网点处高程异常中地形因素的影响ΔNP，无需再计算对应位置布格重力异常中的层间改正与地形改正并用Stokes积分将其转换为高程异常中地形因素的影响。
为了利用先前预留的检核点对所得局部似大地水准面模型进行检验，可以对格网化后的局部似大地水准面模型采用反距离加权插值的方法，计算检核点位置对应的高程异常，并将计算结果与GPS水准数据进行对比，以此评价局部似大地水准面模型的精度。
本实施例中，对本发明提供的顾及重力地形改正的GPS高程拟合方法与常规未经地形改正的GPS高程拟合方法得到的局部似大地水准面模型使用检核点进行检核，将GPS高程拟合方法计算结果与实际GPS水准数据的高程异常相减得到检核点高程异常的拟合残差值，该残差即为拟合结果与真实情况之间的不符合程度，检核点拟合残差的均方根可以视为拟合模型的外符合精度。比较多面函数方法、基于二次曲面的Shepard方法和基于二次曲面的多面函数方法的拟合残差结果，从图4中可以看出，顾及重力地形改正的GPS高程拟合方法的检核点拟合残差相较未经改正的常规GPS高程拟合方法大多有较好的改进，如图4中(a)部分多面函数方法的3号和17号点，图4中(b)部分基于二次曲面的Shepard方法的1号、2号、3号、4号和10号点以及图4中(c)部分基于二次曲面的多面函数方法的1号、3号和17号点，其改进效果更为明显。此外，分别统计三种GPS高程拟合方法未经改正以及经过重力地形改正后检核点高程异常的残差均方根以及极值情况，得到下表：
Figure BDA0003127370690000091
Figure BDA0003127370690000101
从上表可以看出，本发明提供的顾及重力地形改正的GPS高程拟合方法的拟合模型外符合精度相较于未经改正的方法，在多面函数法上有36.15％改进效果，在基于二次曲面的Shepard方法上有52.36％的改进效果，在基于二次曲面的多面函数方法上有50.46％的改进效果。可以认为利用本发明提出的顾及地形改正的GPS拟合方法可以明显改善常规GPS拟合方法在山区的拟合效果。
具体实施时，本发明技术方案提出的方法可由本领域技术人员采用计算机软件技术实现自动运行流程，实现方法的***装置例如存储本发明技术方案相应计算机程序的计算机可读存储介质以及包括运行相应计算机程序的计算机设备，也应当在本发明的保护范围内。
在一些可能的实施例中，提供一种顾及重力地形改正的GPS高程拟合***，包括以下模块，
第一模块，用于根据测区的数字高程模型中包含的大地高信息以及测区格网信息计算格网点处布格重力异常的地形改正TC项和层间改正BP项；
第二模块，用于将布格重力异常的地形改正和层间改正相加，所得结果δg使用Stokes积分公式转换为对似大地水准面的影响，得到格网点处高程异常中地形因素的影响ΔNP；
第三模块，用于选取测区中的GPS水准点，利用测点上的大地高H和正常高h相减高程异常ξ，再去除GPS水准点处地形因素的影响ΔN，得到经过地形改正后高程异常ξ′；
第四模块，用于对地形改正后GPS水准点的高程异常进行数值拟合，得到地形改正后的局部似大地水准面模型N'；
第五模块，用于将第四模块所得地形改正后的局部似大地水准面模型N'与对应位置计算得到的高程异常中的地形影响ΔNP结合，得到真实的局部似大地水准面N，作为最终所得顾及重力地形改正的GPS高程拟合成果。
在一些可能的实施例中，提供一种顾及重力地形改正的GPS高程拟合***，包括处理器和存储器，存储器用于存储程序指令，处理器用于调用存储器中的存储指令执行如上所述的一种顾及重力地形改正的GPS高程拟合方法。
在一些可能的实施例中，提供一种顾及重力地形改正的GPS高程拟合***，包括可读存储介质，所述可读存储介质上存储有计算机程序，所述计算机程序执行时，实现如上所述的一种顾及重力地形改正的GPS高程拟合方法。
本文中所描述的具体实施例仅仅是对本发明精神作举例说明。本发明所属技术领域的技术人员可以对所描述的具体实施例做各种各样的修改或补充或采用类似的方式替代，但并不会偏离本发明的精神或者超越所附权利要求书所定义的范围。

Claims (10)

1.一种顾及重力地形改正的GPS高程拟合方法，其特征在于：用于山区地形的GPS高程拟合，包括步骤如下：
步骤1，根据测区的数字高程模型中包含的大地高信息以及测区格网信息计算格网点处布格重力异常的地形改正TC项和层间改正BP项；
步骤2，将布格重力异常的地形改正和层间改正相加，所得结果δg转换为对似大地水准面的影响，得到格网点处高程异常中地形因素的影响ΔNP；
步骤3，选取测区中的GPS水准点，利用测点上的大地高H和正常高h相减高程异常ξ，再去除GPS水准点处地形因素的影响ΔN，得到经过地形改正后高程异常ξ′；
步骤4，对地形改正后GPS水准点的高程异常进行数值拟合，得到地形改正后的局部似大地水准面模型N'；
步骤5，将步骤4所得地形改正后的局部似大地水准面模型N'与对应位置计算得到的高程异常中的地形影响ΔNP结合，得到真实的局部似大地水准面N，作为最终所得顾及重力地形改正的GPS高程拟合成果。
2.根据权利要求1所述的顾及重力地形改正的GPS高程拟合方法，其特征在于，步骤1中，根据地形起伏与地面点间的关系，直接给出地形改正TC项计算公式如下，
Figure FDA0003689846800000011
其中，G为地心引力常数，ρ(x,y,z)为流动点(x,y,z)的地壳密度，h为流动点的高程，E代表积分区域，r()为距离函数；(xP,yP)为P点平面坐标，hP为P点及其在大地水准面上的投影点P0之间的距离；
层间改正BP项的计算公式为BPP＝-2πGρhP，其中，ρ为地壳密度。
3.根据权利要求2所述的顾及重力地形改正的GPS高程拟合方法，其特征在于：步骤2中，根据重力异常与高程异常间的关系，某地面点P上布格重力异常的层间改正和地形改正对似大地水准面的影响用Stokes积分公式表示成如下形式，
Figure FDA0003689846800000021
其中，γ为平均正常重力，E代表积分区域，R为地球平均半径，
Figure FDA0003689846800000022
为Stokes核函数，
Figure FDA0003689846800000023
和λP为地面点P的纬度和经度，
Figure FDA0003689846800000024
和λ为点P周围流动点的纬度和经度。
4.根据权利要求3所述的顾及重力地形改正的GPS高程拟合方法，其特征在于：步骤2中，布格重力异常的地形改正和层间改正对高程异常间的关系表示为球谐展开式的形式，
Figure FDA0003689846800000025
其中，G为地心引力常数，γ为平均正常重力，R为地球平均半径，ρ为地壳密度，
Figure FDA0003689846800000026
h为流动点
Figure FDA0003689846800000027
的高程；Ynm为n阶m次的完全规格化的球谐函数；
Figure FDA0003689846800000028
和λP为地面点P的纬度和经度，
Figure FDA0003689846800000029
和λ为点P周围流动点的纬度和经度。
5.根据权利要求1所述的顾及重力地形改正的GPS高程拟合方法，其特征在于：步骤3中，根据格网点处高程异常中地形因素的影响
Figure FDA00036898468000000210
以反距离加权插值的方式得到GPS水准点处高程异常中地形因素的影响ΔN，反距离加权插值表示成如下形式，
Figure FDA00036898468000000211
其中，Pi,i＝1,2,3,4为距离该GPS水准点最近的4个格网点，di为GPS水准点至格网点的平面距离。
6.根据权利要求1所述的顾及重力地形改正的GPS高程拟合方法，其特征在于：步骤4中，对地形改正后的GPS水准点高程异常进行数值拟合采用多项式方法、Shepard方法、多面函数方法或BP神经网络方法。
7.根据权利要求1或2或3或4或5或6所述的顾及重力地形改正的GPS高程拟合方法，其特征在于：将步骤4经过数值拟合构建得到的地形改正后的局部似大地水准面模型N'用格网数据形式表达，使得拟合模型的格网点与步骤1所述测区格网一一对应，则在步骤5恢复真实的局部似大地水准面模型N时，直接利用步骤2中的格网点处高程异常中地形因素的影响ΔNP，无需再计算对应位置布格重力异常中的层间改正与地形改正并用Stokes积分转换为高程异常中地形因素的影响。
8.一种顾及重力地形改正的GPS高程拟合***，其特征在于：用于实现如权利要求1-7任一项所述的一种顾及重力地形改正的GPS高程拟合方法。
9.根据权利要求8所述顾及重力地形改正的GPS高程拟合***，其特征在于：包括以下模块，
第一模块，用于根据测区的数字高程模型中包含的大地高信息以及测区格网信息计算格网点处布格重力异常的地形改正TC项和层间改正BP项；
第二模块，用于将布格重力异常的地形改正和层间改正相加，所得结果δg使用Stokes积分公式转换为对似大地水准面的影响，得到格网点处高程异常中地形因素的影响ΔNP；
第三模块，用于选取测区中的GPS水准点，利用测点上的大地高H和正常高h相减高程异常ξ，再去除GPS水准点处地形因素的影响ΔN，得到经过地形改正后高程异常ξ′；
第四模块，用于对地形改正后GPS水准点的高程异常进行数值拟合，得到地形改正后的局部似大地水准面模型N'；
第五模块，用于将第四模块所得地形改正后的局部似大地水准面模型N'与对应位置计算得到的高程异常中的地形影响ΔNP结合，得到真实的局部似大地水准面N，作为最终所得顾及重力地形改正的GPS高程拟合成果。
10.根据权利要求8所述顾及重力地形改正的GPS高程拟合***，其特征在于：包括处理器和存储器，存储器用于存储程序指令，处理器用于调用存储器中的存储指令执行如权利要求1-7任一项所述的一种顾及重力地形改正的GPS高程拟合方法。
CN202110692845.1A 2021-06-22 2021-06-22 一种顾及重力地形改正的gps高程拟合方法及*** Active CN113378471B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110692845.1A CN113378471B (zh) 2021-06-22 2021-06-22 一种顾及重力地形改正的gps高程拟合方法及***

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110692845.1A CN113378471B (zh) 2021-06-22 2021-06-22 一种顾及重力地形改正的gps高程拟合方法及***

Publications (2)

Publication Number Publication Date
CN113378471A CN113378471A (zh) 2021-09-10
CN113378471B true CN113378471B (zh) 2022-09-02

Family

ID=77578386

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110692845.1A Active CN113378471B (zh) 2021-06-22 2021-06-22 一种顾及重力地形改正的gps高程拟合方法及***

Country Status (1)

Country Link
CN (1) CN113378471B (zh)

Families Citing this family (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113985490B (zh) * 2021-09-22 2023-05-05 中国人民解放军战略支援部队信息工程大学 利用地形和地壳密度数据进行地表重力仿真的方法及装置
CN115112090A (zh) * 2022-06-29 2022-09-27 重庆市勘测院((重庆市地图编制中心)) 一种基于少量gnss水准点的山区数字高程基准建立方法
CN117538911A (zh) * 2024-01-10 2024-02-09 武汉大学 一种实时正常高获取方法及***

Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102175213A (zh) * 2011-01-30 2011-09-07 中国测绘科学研究院 Gps水准高程传递方法

Family Cites Families (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US7616449B2 (en) * 2006-09-25 2009-11-10 Appareo Systems, Llc Crash-hardened memory device and method of creating the same
CN111257956A (zh) * 2020-04-02 2020-06-09 吉林省水利水电勘测设计研究院 一种基于Matlab的区域似大地水准面精化方法

Patent Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102175213A (zh) * 2011-01-30 2011-09-07 中国测绘科学研究院 Gps水准高程传递方法

Also Published As

Publication number Publication date
CN113378471A (zh) 2021-09-10

Similar Documents

Publication Publication Date Title
CN113378471B (zh) 一种顾及重力地形改正的gps高程拟合方法及***
CN112556660B (zh) 基于卫星测高数据的海域重力异常反演方法及***
CN111257956A (zh) 一种基于Matlab的区域似大地水准面精化方法
Abdalla et al. The evaluation of the New Zealand's geoid model using the KTH method
CN102175213A (zh) Gps水准高程传递方法
Fryer et al. On the accuracy of heighting from aerial photographs and maps: implications to process modellers
Oluyori et al. Computations of geoid undulation from comparison of GNSS/levelling with EGM 2008 for geodetic applications
Sánchez Strategy to establish a global vertical reference system
KR101076299B1 (ko) 중력 시뮬레이션 자료 구축 방법 및 이를 수행하기 위한 시스템
Nsiah Ababio et al. Compilation of the new detailed geoid model HKGEOID-2022 for the Hong Kong territories
CN108038086A (zh) 一种基于像素尺度的dem数据误差评价与校正方法
CN109446676B (zh) 一种高程***基准面确定方法及***
CN113985490B (zh) 利用地形和地壳密度数据进行地表重力仿真的方法及装置
CN116299440A (zh) 一种基于D-InSAR技术的滑坡防治效果评价方法
Okiemute et al. Relationship between Polynomial Geometric Surfaces Terms and Observation Points Numbers and Effect in the Accuracy of Geometric Geoid Models
Willberg et al. Analysis of GOCE omission error and its contribution to vertical datum offsets in Greece and its Islands
Majkráková et al. An analysis of methods for gravity determination and their utilization for the calculation of geopotential numbers in the Slovak national levelling network
Papp et al. Effect of the difference between surface and terrain models on gravity field related quantities
Kemboi et al. Estimation of orthometric height using EGM2008 and GPS over Nairobi county and its environs
Yilmaz et al. A research of consistencies and progresses of geoid models in Turkey
Mo et al. Evaluation of AW3D30 Elevation Accuracy in China.
Huang et al. An adaptive approach to optimise regional geoid undulation model for engineering applications
Zaletnyik et al. Combination of GPS/leveling and the gravimetric geoid by using the thin plate spline interpolation technique via finite element method
Spiroiu et al. Geo-gravimetric quasi-geoid determination over Romania
Ferrara et al. Adaptation of the Global Geoid Model EGM2008 on Campania Region (Italy) based on geodetic network points

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant