CN102930176B - Nine-grid five-spline interpolation method for geoid model - Google Patents

Nine-grid five-spline interpolation method for geoid model Download PDF

Info

Publication number
CN102930176B
CN102930176B CN201210478131.1A CN201210478131A CN102930176B CN 102930176 B CN102930176 B CN 102930176B CN 201210478131 A CN201210478131 A CN 201210478131A CN 102930176 B CN102930176 B CN 102930176B
Authority
CN
China
Prior art keywords
interpolation
point
node
height anomaly
height
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.)
Expired - Fee Related
Application number
CN201210478131.1A
Other languages
Chinese (zh)
Other versions
CN102930176A (en
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.)
Changsha University of Science and Technology
Original Assignee
Changsha University of Science and Technology
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 Changsha University of Science and Technology filed Critical Changsha University of Science and Technology
Priority to CN201210478131.1A priority Critical patent/CN102930176B/en
Publication of CN102930176A publication Critical patent/CN102930176A/en
Application granted granted Critical
Publication of CN102930176B publication Critical patent/CN102930176B/en
Expired - Fee Related legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Landscapes

  • Position Fixing By Use Of Radio Waves (AREA)
  • Image Generation (AREA)

Abstract

A method utilizing (quasi) geoid grid interpolation for calculating height anomaly of any point belongs to the field of geodesy and measurement engineering technology. GPS (Global Positioning System) geoid height measurement has advantages of high-speed and high-precision, however, the normal height system in which the GPS ground height has to be converted into the normal height is adopted by China. The method provided by the invention, which is based on the (quasi) geoid grid, utilizes the nine-grid five-spline interpolation to calculate the height anomaly. According to the method, based on nine grids and sixteen joints closest to a point for interpolation, four joints on the same straight line are respectively adopted to form polishing spline curves interpolation which are used for calculating the height anomaly of the middle point, and then spline curves formed by the middle points are adopted to calculate the height anomaly of the point for interpolation, and the 1985 Yellow Sea height can be obtained through subtracting the height anomaly by the GPS ground height of the point. The method provided by the invention has the advantages that the continuity and smoothness problems of a current geoid grid interpolation are thoroughly solved, and the interpolation result is sole and is very high in computation accuracy. Therefore, the method can be widely used for modern height measurement.

Description

Geoid's model nine grid five Spline Interpolation Method
Technical field
The present invention is a kind of method calculated for geoid surface or quasigeoid grid model spline interpolation, belongs to " physical geodesy " technical field in " Geodesy and Survey Engineering " subject.
Figure viewed from behind technology
The normal level that the GPS geodetic height that certain any height anomaly is numerically equal to this point deducts this point is high.In Surveying Engineering application, the universal GPS of the making geodetic height of GPS is converted into GPS survey becomes a kind of urgent need of production, and height anomaly is the parameters necessary of this conversion process.The method combined by gravity and geometry, first refines according to tight theory and survey the geoid's model in district, then the height anomaly calculating unknown point is current prefered method.The present invention mainly solves the problem calculating unknown point height anomaly on the basis of geoid surface or quasigeoid grid model.Up to now, this problem has proposed three technical schemes, but these three technical schemes all exist deficiency in various degree.As everyone knows, geoid surface is a continuous print, smooth, complicated curved surface.First technical scheme carries out unknown point height anomaly interpolation with inverse distance-weighting in single geoid surface graticule mesh, and this method exists serious problem, and namely interpolation curved surface gross distortion can occur, and curved surface disconnects, discontinuous, rough; Second technical scheme is a kind of easy interpolation method, during its interpolation calculation remains and carries out in single graticule mesh, be not but adopt inverse distance weighted interpolation, but adopt linear interpolation, this method also has problems, namely interpolation result is not smooth curved surface, but pinch plane, there is obvious cut-off rule in mesh edge, is similar to folding line, and cusp is formed at node place, cause curved surface rough; 3rd technical scheme is put forward by applicant, the interpolation calculation of the method is carried out in four graticule mesh, adopt Biquadratic Polynomial matching, the height anomaly of unknown point is obtained through continuous cubic interpolation, but in use, we find that Biquadratic Polynomial matching exists " imperial lattice phenomenon ", and namely curve exists unnecessary swing, and interpolation result may not uniquely, cause reducing Security and feasibility in use.Substantially be all adopt these three kinds of technical schemes at present both at home and abroad.
Patent of the present invention is the 4th technical scheme that we propose, its interpolation is carried out in nine graticule mesh, interpolation curve is five smoothing splint curves, its advantage is the continuity and the slickness problem that solve geoid surface graticule mesh interpolation up hill and dale, interpolation result is unique, and computational accuracy is very high.
Summary of the invention
A kind of geoid's model nine grid five Spline Interpolation Method, is characterized in that, comprise the following steps:
L) regularization of geoid surface grid model form:
(1.1) regularization of indication of the present invention, comprise grid model from east to west, by north orientation south equidistant arrangement, graticule mesh is neat continuously;
(1.2) lower left corner of graticule mesh, the terrestrial coordinate in the upper right corner is indicated;
(1.3) spacing of graticule mesh east-west direction and the spacing of North and South direction is indicated;
2) by interpolation point terrestrial coordinate, the graticule mesh participating in calculating is determined:
(2.1) according to the coordinate of interpolation point, determine which graticule mesh interpolation point drops into, and writes down its numbering i;
(2.2) eight graticule mesh adjacent with graticule mesh i are determined;
(2.3) write down this step 2) in (2.1), (2.2) determined nine graticule mesh node numbering;
3) coordinate, the height anomaly of 16 numbering nodes of nine graticule mesh is determined; Wherein, the numbering of 16 nodes of nine graticule mesh from top to bottom, from left to right number consecutively be 1,2,3 ..., 16;
4) by its height anomaly of interpolation point Geodetic Coordinate Calculation, concrete steps are as follows:
(4.1) cross the coordinate of interpolation point, parallel parallel draws parallel lines;
(4.2) four nodes a, b, c, d that parallel lines are crossing with graticule mesh are determined;
(4.3) by the grid line at a, b, c, d tetra-node places, spline curve fitting is formed;
(4.3.1) by equidistant 4 structure minimum curvature battens, parameters, comprises the number of times of splines number of times, batten correction; Read in the coordinate of known node, and be " batten number of times adds 1 " individual known point at end points place to extension number, the coordinate of extension point is determined according to being calculated by two or a consecutive point coordinate;
(4.3.2) set 4 nodes of wish formation curve as r i(x i, y i), i=1 ~ 4, then the δ spline function curve connecting k batten k-1 the polishing of each node is: wherein r ifor the height anomaly of node; R kh () is the height anomaly of some h matching polishing, put h and 4 node on the same line, and 1≤h≤4;
(4.3.3) Ω kx () is δ Gene radiotherapy: Ω k ( x ) = Σ j = 0 k + 1 ( - 1 ) j k + 1 j ( x + k + 1 2 - j ) + k 1 k ! , Character expression wherein calculates by following: n m = n ! ( n - m ) ! m ! , Z + k = Z k , Z > 0 0 , Z ≤ 0
(4.3.4) height anomaly of polishing match point is calculated, and matching is obtained comparing at the calculation level height anomaly of site position and the original height anomaly of node, if less than original node abnormal elevation value, then its difference is added on original height anomaly; Otherwise, difference is reduced from the original height anomaly of node, namely the original height anomaly of node is revised;
(4.3.5) again calculate final polishing match point height anomaly with revised original node height anomaly, obtain than the front polishing matched curve once more approaching original node;
(4.4) in nine graticule mesh, by 1,2,3,4 four node structure Article 1 SPL, the height anomaly of interpolation a point; By 5,6,7,8 four node structure Article 2 SPL, the height anomaly of interpolation b point; By 9,10,11,12 4 node structure Article 3 SPL, the height anomaly of interpolation c point; By 13,14,15,16 4 node structure Article 4 SPL, the height anomaly of interpolation d point;
(4.5) according to the height anomaly of a, b, c, d tetra-nodes, according to (4.3.1) in step (4.3) to the step in (4.3.5), structure Article 5 SPL, matching a, b, c, d tetra-nodes, the height anomaly of interpolation P point.
2., because geoid surface grid model and quasigeoid grid model have identical form, its difference is only node height anomaly slightly difference.Therefore, geoid's model nine grid five Spline Interpolation Method according to claim 1, also can be applicable in quasigeoid grid model interpolation calculation.
Accompanying drawing explanation
Fig. 1 geoid's model nine grid five Spline Interpolation Method schematic diagram
Embodiment
Geoid's model nine grid five Spline Interpolation Method, its elementary object is based on grid model, is gone out the abnormal elevation value of position, arbitrfary point by the abnormal elevation value interpolation calculation of node.Two innovative points of the present invention are: one is that the world proposes employing nine graticule mesh first to carry out interpolation calculation at home, completely avoid the unique and rough discontinuous problem of curved surface of interpolation result; Two is in the calculating of a unknown point height anomaly, and introduced by five smoothing splint curves in the interpolation calculation of geoid surface grid model, interpolation result is reliable, and interpolation precision is very high.Principle constructs smooth SPL according to 4 known points equidistantly distributed, and obtained the height anomaly of impact point position by five interpolation.Comprise the following steps in concrete enforcement:
One, the regularization of geoid surface grid model form;
The grid model regularization of indication of the present invention, refers to that grid model forms the file layout of specifying: the lower left corner of graticule mesh, the terrestrial coordinate in the upper right corner and graticule mesh east-west direction spacing and North and South direction spacing are indicated in file first trip; Second row starts the abnormal elevation value for node, and often row from east to west, often arranges by the north orientation equidistant arrangement in south, and graticule mesh is neat continuously.
Two, by interpolation point terrestrial coordinate, the height anomaly of interpolation point is calculated;
1) if interpolation point terrestrial coordinate overlaps with certain node, then its height anomaly is exactly the height anomaly of node.
2) if interpolation point falls on the line of certain two node, then the SPL that its height anomaly is formed by four nodes is determined, these four nodes are located along the same line, and nearest with interpolation point.
3) if interpolation point falls into graticule mesh blank space, then determine the graticule mesh of participating in calculating, according to the coordinate of interpolation point, determine which graticule mesh interpolation point drops into, and writes down its numbering i; Determine eight graticule mesh adjacent with graticule mesh i; Write down the numbering of nine mesh node; Determine coordinate and the height anomaly of 16 numbering nodes of nine graticule mesh.
4) interpolation calculation height anomaly according to the following steps, implementation step is as follows:
(4.1) cross the coordinate of interpolation point, parallel parallel draws parallel lines, and these parallel lines stride across three nearest graticule mesh of interpolation point, and intersect at four intersection points a, b, c, d with graticule mesh;
(4.2) four collinear and nearest with putting a nodes form spline fit curve, thus determine the height anomaly of an a; In like manner obtain the height anomaly of other three intersection points such as some b, c, d etc.
5) spline curve construction step is as follows:
(5.1) by equidistant 4 structure minimum curvature battens, parameters, comprises the number of times of splines number of times, batten correction; Read in the coordinate of known node, and be " batten number of times adds 1 " individual known point at end points place to extension number, the coordinate of extension point calculates according to two consecutive point coordinates to be determined;
(5.2) set the node of wish formation curve as r i(x i, y i), i=1 ~ 4, then the δ spline function curve connecting k batten k-1 the polishing of each node is: wherein r ifor the height anomaly of node; R kh () is the height anomaly of some h matching polishing, put h and 4 node on the same line, and 1≤h≤4;
(5.3) Ω kx () is δ Gene radiotherapy: Ω k ( x ) = Σ j = 0 k + 1 ( - 1 ) j k + 1 j ( x + k + 1 2 - j ) + k 1 k ! , Character expression wherein calculates by following: n m = n ! ( n - m ) ! m ! , Z + k = Z k , Z > 0 0 , Z ≤ 0
(5.4) height anomaly of polishing match point is calculated, and matching is obtained comparing at the calculation level height anomaly of site position and the original height anomaly of node, if less than original node abnormal elevation value, then its difference is added on original height anomaly; Otherwise, difference is reduced from the original height anomaly of node, namely the original height anomaly of node is revised;
(5.5) again calculate final polishing match point height anomaly with revised original node height anomaly, obtain than the front polishing matched curve once more approaching original node.
6) as shown in Figure of description, in nine graticule mesh, by 1,2,3,4 four node structure Article 1 SPL, the height anomaly of interpolation a point; By 5,6,7,8 four node structure Article 2 SPL, the height anomaly of interpolation b point; By 9,10,11,12 4 node structure Article 3 SPL, the height anomaly of interpolation c point; By 13,14,15,16 4 node structure Article 4 SPL, the height anomaly of interpolation d point.Interpolation point P and a, b, c, d tetra-intersection points are on the same line, according to the height anomaly of a, b, c, d tetra-nodes, according to step 5) in the method that describes to (5.5) of sub-step (5.1), structure Article 5 SPL, the height anomaly of interpolation P point.

Claims (2)

1. geoid's model nine grid five Spline Interpolation Method, is characterized in that, comprise the following steps:
L) regularization of geoid surface grid model form:
(A) regularization of indication of the present invention, comprise grid model from east to west, by north orientation south equidistant arrangement;
(B) lower left corner of graticule mesh, upper right corner terrestrial coordinate is indicated;
(C) each spacing of graticule mesh east-west direction and the spacing of North and South direction is indicated;
2) by interpolation point terrestrial coordinate, the graticule mesh of participating in calculating is determined:
I (), according to the coordinate of interpolation point, determines which graticule mesh interpolation point drops into, and writes down its numbering i;
(ii) eight graticule mesh adjacent with graticule mesh i are determined;
(iii) write down this step 2) in (i), (ii) nine graticule mesh the numbering of node;
3) coordinate, the height anomaly of 16 nodes of nine graticule mesh is determined; Wherein, the numbering of 16 nodes of nine graticule mesh from top to bottom, from left to right number consecutively be 1,2,3 ..., 16;
4) by its height anomaly of interpolation point Geodetic Coordinate Calculation, concrete steps are as follows:
(A) cross the coordinate of interpolation point, parallel parallel draws parallel lines;
(B) four nodes a, b, c, d that parallel lines are crossing with graticule mesh are determined;
(C) by the grid line at a, b, c, d tetra-node places, by following (C-1) to (C-5) five steps, spline curve fitting is formed;
(C-1) by equidistant 4 structure minimum curvature battens, parameters, comprises the number of times of splines number of times, batten correction; Read in the coordinate of known node, and be " batten number of times adds 1 " individual known point at end points place to extension number, the coordinate of extension point calculates according to the coordinate of two consecutive point to be determined;
(C-2) set 4 nodes of wish formation curve as r i(x i, y i), i=1 ~ 4, then the δ spline function curve connecting k batten k-1 the polishing of each node is: wherein r ifor the height anomaly of node; R kh () is the height anomaly of some h matching polishing, put h and 4 node on the same line, and 1≤h≤4;
(C-3) Ω kx () is k δ Gene radiotherapy: Ω k ( x ) = Σ j = 0 k + 1 ( - 1 ) j k + 1 j ( x + k + 1 2 - j ) + k 1 k ! , Character expression wherein calculates by following: n m = n ! ( n - m ) ! m ! , Z + k = Z k , Z > 0 0 , Z ≤ 0 ;
(C-4) height anomaly of polishing match point is calculated, and matching is obtained comparing at the calculation level height anomaly of site position and the original height anomaly of node, if less than original node abnormal elevation value, then its difference is added on original height anomaly; Otherwise, difference is reduced from the original height anomaly of node, namely the original height anomaly of node is revised;
(C-5) with the polishing match point height anomaly that revised original node height anomaly iterative computation is final, obtain than the front polishing matched curve once more approaching original node;
(D) in nine graticule mesh, by 1,2,3,4 four node structure Article 1 SPL, the height anomaly of interpolation a point; By 5,6,7,8 four node structure Article 2 SPL, the height anomaly of interpolation b point; By 9,10,11,12 4 node structure Article 3 SPL, the height anomaly of interpolation c point; By 13,14,15,16 4 node structure Article 4 SPL, the height anomaly of interpolation d point;
(E) according to the height anomaly of a, b, c, d tetra-nodes, according to (C-1) in step (C) to the step in (C-5), structure Article 5 SPL, the height anomaly of interpolation P point.
2. geoid's model nine grid five Spline Interpolation Method according to claim 1, it is characterized in that, the method can be applicable in quasigeoid grid model interpolation calculation.
CN201210478131.1A 2012-11-22 2012-11-22 Nine-grid five-spline interpolation method for geoid model Expired - Fee Related CN102930176B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201210478131.1A CN102930176B (en) 2012-11-22 2012-11-22 Nine-grid five-spline interpolation method for geoid model

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201210478131.1A CN102930176B (en) 2012-11-22 2012-11-22 Nine-grid five-spline interpolation method for geoid model

Publications (2)

Publication Number Publication Date
CN102930176A CN102930176A (en) 2013-02-13
CN102930176B true CN102930176B (en) 2015-03-25

Family

ID=47644973

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201210478131.1A Expired - Fee Related CN102930176B (en) 2012-11-22 2012-11-22 Nine-grid five-spline interpolation method for geoid model

Country Status (1)

Country Link
CN (1) CN102930176B (en)

Families Citing this family (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104933291B (en) * 2015-05-06 2017-08-25 中国石油大学(华东) Altimetry Data mean recovery time production method based on net―function
US10070259B1 (en) * 2017-02-24 2018-09-04 Here Global B.V. Altitude map for indoor positioning services
CN111366130A (en) * 2020-03-03 2020-07-03 山东交通学院 Elevation anomaly interpolation and estimation method

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101319893A (en) * 2008-07-18 2008-12-10 东南大学 Method for accurately determining region height anomaly
CN101839710A (en) * 2010-06-12 2010-09-22 中国测绘科学研究院 Method for optimizing quasi-geoid calculation
CN102607513A (en) * 2012-02-14 2012-07-25 东南大学 Method for carrying out quasigeoid refining on superlarge region on basis of seamless partitioning technology
WO2012138837A1 (en) * 2011-04-08 2012-10-11 Fleetmatics Irl Limited System and method for providing an electronic representation of a route

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101319893A (en) * 2008-07-18 2008-12-10 东南大学 Method for accurately determining region height anomaly
CN101839710A (en) * 2010-06-12 2010-09-22 中国测绘科学研究院 Method for optimizing quasi-geoid calculation
WO2012138837A1 (en) * 2011-04-08 2012-10-11 Fleetmatics Irl Limited System and method for providing an electronic representation of a route
CN102607513A (en) * 2012-02-14 2012-07-25 东南大学 Method for carrying out quasigeoid refining on superlarge region on basis of seamless partitioning technology

Non-Patent Citations (4)

* Cited by examiner, † Cited by third party
Title
On the non-uniqueness of local quasi-geoids computed from terrestrial gravity anomalies;I.PRUTKIN,ET AL.;《J GEOD》;20081231;第82卷(第3期);147~156 *
似大地水准面格网双二次多项式插值方法;邓兴升等;《测绘学报》;20090228;第38卷(第1期);35~40 *
厘米级高程异常地形影响的算法及特征分析;章传银等;《测绘学报》;20061130;第35卷(第4期);308~314 *
移动格林基函数样条二维插值算法研究;邓兴升等;《大地测量与地球动力学》;20111231;第31卷(第6期);69~72 *

Also Published As

Publication number Publication date
CN102930176A (en) 2013-02-13

Similar Documents

Publication Publication Date Title
CN105843234B (en) A kind of two-dimentional Route planner that UUV detours to round barrier geometry
CN103927788A (en) Building ground feature DEM manufacturing method based on city vertical planning
CN105138788A (en) Full-parameterization establishing method for road three-dimensional model
CN103217688B (en) Airborne laser radar point cloud adjustment computing method based on triangular irregular network
CN103927418A (en) Method for manufacturing meshed drainage channels of urban road drains on basis of DEM (digital elevation model)
CN105355042B (en) A kind of road network extraction method based on taxi GPS
CN102879820B (en) Based on the three-dimensional table layer model construction method of triangle gridding
CN102930176B (en) Nine-grid five-spline interpolation method for geoid model
CN105180811A (en) Laser scanner calibration method, based on ground objects with characteristics of the same name, for mobile measuring system
CN107918957A (en) A kind of three-dimensional building model simplifying method of holding structure and textural characteristics
CN101838958A (en) Road gradient detection method
CN103400416A (en) City environment robot navigation method based on multi-layer probabilistic terrain
CN106640092A (en) Posture rechecking method for shield machine used for tunnel construction
CN111005273A (en) Temporary road arrangement method for construction
CN103970837A (en) Discontinuous DEM classified manufacturing method based on urban land and vertical planning
CN103455709B (en) A kind of super-resolution method for digital elevation model and system thereof
CN109297501A (en) A kind of method that GIS single line road network automatically generates two-wire road network
CN106569260B (en) For high-dip structure area two-dimension earthquake speed change sky school drawing methods
CN104700710A (en) Simulation map for house property mapping
CN109166173B (en) Multi-precision three-dimensional mapping data fusion method based on BIM
CN104123671B (en) Meteorological reference station partition method
CN105466391A (en) A tower base section generating method by utilization of a digital elevation model and field-data correction
CN105224748A (en) A kind of section preprocess method of non-uniform beam finite element model
CN108596397A (en) One kind being based on hydrological model multimachine flight course planning method
CN105203070A (en) Calculation method for soil deformation caused by double-circle shield tunnel construction deflection

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
C14 Grant of patent or utility model
GR01 Patent grant
CF01 Termination of patent right due to non-payment of annual fee
CF01 Termination of patent right due to non-payment of annual fee

Granted publication date: 20150325

Termination date: 20181122