GNSS Positioning by CORS and EGM2008 in Jilin Province, China

The Continuously Operating Reference Station (CORS) technique has been widely applied in land resource management, surveying, mapping, deformation monitoring, precise navigation, etc. This article analyzed the positioning method using EGM2008 and CORS of Jilin Province, China. The vertical transformation of EGM2008 from WGS84 to China’s CGCS2000 datum and the horizontal coordinate transformation from CGCS2000 to a triangulation coordinate system were discussed. The results indicated that a local geoid with respect to CGCS2000 can be transferred from EGM2008 with the same accuracy, and the geoid correction between CGCS2000 and WGS84 varied from 0.023 m to 0.111 m. The coordinate transformation method based on the curve surface approximation method indicated that the theoretical error was less than 0.09 m in the grid within 10° longitudinal and 5° latitudinal, and less than 0.3 m in large area and 0.1 m in small area in field validation. The method proposed in this article expanded the positioning result and its application for JLCORS and other CORS with local datum.

reference stations with similar functions as the national geographic information benchmark, and provide services including GPS differential positioning, navigation, deformation, etc. [7]. In Asia, The Geographical Survey Institute of Japan (GSI) established in 1993 the GPS observational network, which possesses more than 1200 continuous GPS stations with a mean separation distance of 25 km covering the whole country, and can be applied in crustal deformation monitoring, prediction of earthquakes, meteorology, GPS real-time positioning, etc. [8].
China's first CORS was the Shenzhen Continuous Operational Reference System (SZCORS) established in 2003 by the local government of Shenzhen, in Guangdong Province. SZCORS includes five GNSS continuously operating reference stations and a local geoid of 1 cm accuracy, and provides positioning services (3 cm on horizontal and 5 cm on vertical) in support of geodesy, engineering surveying, meteorological monitoring, earthquake monitoring, precision navigation, etc. [9]. Since then Beijing, Tianjin, Shanghai, Chengdu, Chongqing, Kunming, Wuhan, Hong Kong and some other cities also established CORS [10].
The Earth Gravitational Model 2008 (EGM2008) was released by the EGM Development Team of the National Geospatial-Intelligence Agency (NGA). It is a high order gravitational model complete to spherical harmonic degree and order 2159, while coefficients extending to degree 2190 and order 2159 was also included [11][12][13]. The overall accuracy of EGM2008 in China is about 20-24 cm in Western China, 12 cm in Mid-east China, and 9 cm in North China [14,15].
Jilin Province is located in Northeast China. The terrain inclines from southeast to northwest, with mountains in the eastern part and plains in the middle and western part. The average elevation is 260 m (Figure 1). The Jilin Continuously Operating Reference Station (JLCORS) was established by the Bureau of Surveying Mapping and Geoinformation of Jilin Province. JLCORS consists of 49 reference stations covering the entire Jilin Province (121˝37 1 -131˝18 1 E, 40˝51 1 -46˝19 1 N, 187,400 square kilometers), and the present service aims to provide high accuracy (2-3 cm horizontal and 4-5 cm vertical) coordinates with CGCS2000 datum [16][17][18]. Virtual Reference Station (VRS) technology was applied in JLCORS, and consists of four major sections including a continuous GPS reference station network, data processing center, data communication section and user section [19], which integrates internet technology, wireless communication technology, computer networks and GPS positioning technology. JLCORS provides coordinates within CGCS2000, while many users need the coordinates with respect to old coordinate systems based on triangulation technology such as Beijing54 or Xi'an80 datum [20][21][22], which means the conversion problem from CGCS2000 to these two old datum formats needs to be solved. Another issue that needs to be solved for the application of JLCORS is the definition of local geoid, which is independent in the height conversion from ellipsoid to normal height or orthometric height (the difference of these two height systems is too small to discriminate considering the 200 m average height of the study area). This article studied the coordinate transformation model from geocentric CGCS2000 to the old reference ellipsoid systems (Beijing54 or Xi'an80), and calculated the local geoid based on EGM2008. The accuracy of coordinate transformation and the local geoid was analyzed using field data, and the application of the proposed method is discussed. orthometric height (the difference of these two height systems is too small to discriminate considering the 200 m average height of the study area). This article studied the coordinate transformation model from geocentric CGCS2000 to the old reference ellipsoid systems (Beijing54 or Xi'an80), and calculated the local geoid based on EGM2008. The accuracy of coordinate transformation and the local geoid was analyzed using field data, and the application of the proposed method is discussed.

Local Geoid Calculation Based on EGM2008
The Earth gravity model provides the information to transfer the ellipsoid height to orthometric height or normal height [23,24]. The spherical harmonic function expansions of geoid undulation can be expressed as [25][26][27][28][29]:

Local Geoid Calculation Based on EGM2008
The Earth gravity model provides the information to transfer the ellipsoid height to orthometric height or normal height [23,24]. The spherical harmonic function expansions of geoid undulation can be expressed as [25][26][27][28][29]: N " R n max ÿ n"2 n ÿ m"0`C nm Cos mλ`S nm Sin mλ˘P nm pSin ϕq (1) where n max is the highest order of the expansion, C nm and S nm are the full normalization coefficients, λ is the geodetic longitude, ϕ is geodetic latitude, and R is the Earth's average radius. Two kinds of bias need to be corrected to convert the EGM2008 geoid from WGS84 to CGCS2000. One is the constant gravity datum difference between China and EGM2008, the other is the ellipsoid difference between WGS84 and CGCS2000 ( Figure 2).
where is the highest order of the expansion, ̅ and ̅ are the full normalization coefficients, is the geodetic longitude, is geodetic latitude, and R is the Earth's average radius. Two kinds of bias need to be corrected to convert the EGM2008 geoid from WGS84 to CGCS2000. One is the constant gravity datum difference between China and EGM2008, the other is the ellipsoid difference between WGS84 and CGCS2000 ( Figure 2). The constant bias of vertical datum between China and EGM2008 can be expressed as: where 84 is the ellipsoid height of WGS84, ℎ is the orthometric height, 84 is the geoid undulation of EGM2008 with respect to WGS84. The ellipsoid datum difference between WGS84 and CGCS2000 can be expressed as: The constant bias of vertical datum between China and EGM2008 can be expressed as: where H 84 is the ellipsoid height of WGS84, h is the orthometric height, N 84 is the geoid undulation of EGM2008 with respect to WGS84. The ellipsoid datum difference between WGS84 and CGCS2000 can be expressed as: where H 2000 is the ellipsoid height with CGCS2000. The corrected geoid undulation within CGCS2000 is expressed as follows: The correction method for h is calculated through the observation of precise levelling points measured by static GPS observation within ITRF2008, and this constant is 0.293 min for Jilin Province [27]. To get the correction of the systematic bias from different coordinate systems, a Trimble GNSS R8 receiver was used by synchronous static observation (8 h) in Changchun, Tonghua, Baicheng and Yanji, which nearly covered the entire province ( Figure 1). The coordinates with CGCS2000 of these points was also measured per hour. The static observation data obtained from these four GPS receivers was calculated in ITRF2008. Transformation parameters from ITRF2008 to CGCS2000 were calculated by the Bursa model. Coordinates of grid points at the surface of WGS84 (zero ellipsoid) was transferred to coordinates in CGCS2000 by a 7-parameter datum transformation model to get the ellipsoid difference between the two datum in CGCS2000: where H 2000 is the ellipsoid height within CGCS2000 transferred from ellipsoid surface with respect to WGS84. The errors of geoid due to the value difference of coordinate with respect to the two geodetic datum can be neglected because this difference is very slight ( Table 1).

Coordinate Transformation
The 7-parameter datum transformation model was applied to transfer coordinate from CGCS2000 to Beijing54 or Xi'an80: where X α , Y α , Z α , X β , Y β and Z β represent the geodetic rectangular coordinates of corresponding points in different coordinate systems, while X, Y, Z, ε X , ε Y , ε Z and k are the seven parameters used in the transformation between different coordinate systems [30][31][32].
This model needs to input the geodetic coordinates of both coordinate systems, while the Beijing54 or Xi'an80 datum was established by triangulation technology, so it is impossible to get a high precision ellipsoid height, which means the rectangular space coordinate has a low accuracy, and it is not feasible for a relatively high accuracy transformation [33,34]. The target for the application of this model is to transfer the horizontal coordinates from CGCS2000 to the local ellipsoid, so in order to apply the model for horizontal coordinate transformation, we modified the input data and the algorithm of the Bursa model based on the curve surface approximation. The CGCS2000 ellipsoid and local ellipsoid are geometrical similar spatial curved surfaces for a relatively small area, and the transformation between the two datum in this small area can be defined as the spatial curve surface approximation process from one ellipsoid surface to another ellipsoid surface. Because the form of both ellipsoids was defined, three common points with the same horizontal coordinates in both of the datum are the minimum point amount needed for the approximation. The algorithm is to set the ellipsoid height as zero for both the CGCS2000 and triangulation local datum firstly, and the equation was solved by an iteration method. The process replaces the geodetic height both in CGCS2000 (H α ) and Beijing54 (H β 0 ) with 0 before the parameter solution, next, it translates Hα to the geodetic height in Beijing54 (denoted by H β 1 ) with the use of the seven parameters above, let d 1 " ř´H α´H 1 β¯2 , and then, use H 1 β to solve the seven parameters again to get the new geodetic height in Beijing54 (denoted by H 2 β ), let d 2 " ř´H α´H 2 β¯2 , and the iteration should continue until d 2 ě d 1 . This method does not have the general spatial seven parameters significant for datum transformation, they are the parameters for the approximation of the ellipsoid surface in the study area. In fact, the ellipsoid height of common points are different in the two coordinate systems, and this difference value was calculated and used to measure the approximation level of the two ellipsoid surfaces during the iteration process. This algorithm does not change the shape or the curvature of the ellipsoid surface, which means that the common points in two ellipsoids cannot totally overlap and there are errors in the common points, and these errors can be attributed as the model error. The model errors change with different ellipsoid parameters or study areas, and can be calculated based on the ellipsoid parameters ( Table 2). The curved surface approximation make it possible to calculate transformation parameters in a situation lacking Beijing54 geodetic height data in the application of the Bursa model. In a large area (longitudinal span is 9˝, and latitudinal span is 4˝in middle latitude zones), the theoretical model error is less than 9 cm. Errors in the common points lead to the oscillation of ellipsoid height correction values in small range. The iteration stopped when d 1 < d 2 , that is, the stop condition of the iteration is the beginning of the oscillation. The approaching of the two ellipsoid surface in the calculation process can be attributed as the distance minimization of common points on the two ellipsoid surfaces. When the iteration stopped, the two ellipsoid surfaces can be attributed as approached by the distance minimization of common points and these distances approach the average vertical distance between the two ellipsoid surfaces at these common points.

The Calculation of Transformation Parameters
We applied 7-parameter datum transformation model to calculate the seven parameters which are used to calculate ∆N. Three field measurement points (Tonghua, Baicheng and Yanji, Figure 1) were set as control points, Changchun (Figure 1) was set as checkpoint during the transformation between CGCS2000 coordinates and WGS84 coordinates. The conversion accuracy of each point is shown in Table 3. The data showed a good compatibility between CGCS2000 and ITRF2008. Statistics of the eight measured values (measurement interval is one hour) in every point by JLCORS shows that the maximum error in x, y was 0.016 m, 0.007 m, respectively, and corresponding RMS was 0.006 m, 0.003 m, which proved the stability of JLCORS.
The systematic biases (∆N) between WGS84 and CGCS2000 obtained by the coordinate transformation increased from southeast to northwest in Jilin Province. The minimum value is 0.023 m, the maximum value is 0.111 m and the average is 0.068 m (Figure 3), and the Geoid undulation of Jilin Province with respect to CGCS2000 is shown in Figure 4.
Sensors 2015, 15 7 CGCS2000 coordinates and WGS84 coordinates. The conversion accuracy of each point is shown in Table 3. The data showed a good compatibility between CGCS2000 and ITRF2008. Statistics of the eight measured values (measurement interval is one hour) in every point by JLCORS shows that the maximum error in x, y was 0.016 m, 0.007 m, respectively, and corresponding RMS was 0.006 m, 0.003 m, which proved the stability of JLCORS.
The systematic biases (∆N) between WGS84 and CGCS2000 obtained by the coordinate transformation increased from southeast to northwest in Jilin Province. The minimum value is 0.023 m, the maximum value is 0.111 m and the average is 0.068 m (Figure 3), and the Geoid undulation of Jilin Province with respect to CGCS2000 is shown in Figure 4.

The Validation of Geoid of Jilin Province with Respect to CGCS2000
Four uniformly distributed GPS leveling points in the plain area and one GPS leveling point in a mountain area (Figure 1) of Jilin Province was used to validate the calculated local geoid, and the result is listed in Table 4. The general accuracy is less than 0.030 m in the plain area and 0.050 m in the mountain area.

Coordinate Transformation from CGCS2000 to Local Geodetic Coordinate
Using the curved surface approximation method, two scales of areas in Jilin Province were selected to calculate the transformation parameters from CGCS2000 to Beijing54, and the resulting transformation accuracy is shown in Table 5. The results indicate that for a large area (50,000 square kilometers), the horizontal errors are less than 0.300 m, while for small areas (250 square kilometers), the horizontal errors are less than 0.100 m, which indicates that the proposed method has fine suitability for the coordinate transformation.

The Validation of Geoid of Jilin Province with Respect to CGCS2000
Four uniformly distributed GPS leveling points in the plain area and one GPS leveling point in a mountain area (Figure 1) of Jilin Province was used to validate the calculated local geoid, and the result is listed in Table 4. The general accuracy is less than 0.030 m in the plain area and 0.050 m in the mountain area.

Coordinate Transformation from CGCS2000 to Local Geodetic Coordinate
Using the curved surface approximation method, two scales of areas in Jilin Province were selected to calculate the transformation parameters from CGCS2000 to Beijing54, and the resulting transformation accuracy is shown in Table 5. The results indicate that for a large area (50,000 square kilometers), the horizontal errors are less than 0.300 m, while for small areas (250 square kilometers), the horizontal errors are less than 0.100 m, which indicates that the proposed method has fine suitability for the coordinate transformation.

Discussion and Conclusions
The Earth Gravitational Model 2008 is a geoid model that has a good adaptability in the plain areas of Northeast China, with high accuracy on a global range. The maximum error in CGCS2000 was 0.030 m and the mean square error was 0.022 m after the system deviation and constant deviation were added to the model. The results showed that it is an efficient and powerful method for improving the accuracy and reliability of GPS height measurements by JLCORS.
The results of applying the transformation algorithm between the CGCS2000 and Beijing54 coordinates showed that the use of the surface approximation method to calculate the seven two-dimensional parameters is a highly effective and simple away of obtaining objective coordinates under conditions of unknown or low accuracy geodetic elevation. Transformation parameters should be calculated by partition if accuracy is highly expected for a large area.
WGS84 coordinates have good compatibility with the CGCS2000 coordinate system adopted in JLCORS that runs 24 h a day without break. With better observation conditions and in the signal coverage, rover stations can achieve high location accuracy utilizing one GPS receiver. In addition, coordinates of object points in different coordinate systems can be rapidly obtained by coordinate transformation, and orthometric height or normal height can also obtained quickly with the application of local geoids. With their characteristics of high efficiency, high accuracy, high reliability and low cost, CORS will play a more significant role in the field of control survey, cadastral survey, land management, water resource management, urban planning, line planning, field data collection, topographic mapping, point laying-out, etc.