Introduction

In recent years, China's BeiDou navigation satellite system (BDS) has made great progress and is expected to achieve the goal of global coverage by 2020. The coordinate system used is the BeiDou coordinate system (BDCS) is a coordinate system which is defined in the framework of the international terrestrial reference frame (ITRS) as the current world geodetic coordinate system (WGS84) (^{Li & Huang, 2016}). Therefore, the research on elevation fitting of global positioning system (GPS) is helpful for the effective conversion between Beidou height and normal height in the future.

At present, through the research of some scholars, the methods of fitting GPS elevation data can be roughly divided into two categories. One is the traditional function fitting method, and the other is the new model fitting method. The traditional methods include weighted average method, plane fitting method, quadratic polynomial surface fitting method, moving quadratic polynomial surface fitting method and cubic polynomial surface fitting method, and the new models include neural network fitting method, support vector machine fitting method and curve cloud model fitting method (^{Zhang et al., 2015}; ^{Ren et al., 2015}; ^{Guo, 2018}; ^{Tang et al., 2016}). In the actual operation, the new fitting method is not suitable for large-scale popularization because of its complex model, long data processing cycle and high requirements for hardware and software. However, the traditional fitting method only considers the relationship between the elevation anomaly and the plane coordinates x and y of the point, and the fitting accuracy is low, even in the practical engineering application, it can not achieve the experimental accuracy. In this paper, based on the traditional GPS elevation fitting method, the geodetic height is included in the solution system, and a GPS elevation fitting model based on ternary polynomial regression is proposed (^{Zhou et al., 2016}; Ding & Sun, 2016; ^{Liu et al., 2016}; ^{Yang & Zie, 2017}; Sun et al., 2014).

Interpretation of elevation anomaly and geodetic height

In Figure 1, the point P and the point P' are on the same vertical line and have the same plane position. Then there are

Figure 1 Elevation system diagram 

In the formulas, H _g is the positive height I _pQ of point P, N is the geoid undulation l _QV , H is the geodetic height I _PU , H _g is the positive height I _P'Q of point P', N' is the geoid undulation I _QV , H' is the geodetic height Ip'w. It can be seen from Figure 1 that N = N', so there is

That is, the positive elevation difference between the P point and the P' point is equal to the geodetic height difference. In the formula: ∇H _g is the positive elevation difference between point P and point P', and ∇H is the geodetic height difference between point P and point P'.

According to the literature (figure 1), the difference between the elevation anomaly η and the geoid N meets

In the formula, g _m is the real average gravity value on the vertical line between the geoid and the earth's surface, and y _m is the average gravity value from the reference ellipsoid along the normal line to the approximate earth's surface (^{Li & Huang, 2016}). For point P’, there is

For formula (4) and (5), there are H _g ≠ H’ _g ,g _m ≠ g’ _m and y _m ≠ y' _m . According to the complexity of the leveling surface and the earth's gravity field, and combined with the formula N = N’,

η ≠ η' can be obtained (^{Shijun & Li, 2015}; ^{Gong et al., 2014}; Škrekovski et al., 2019; ^{Zhao et al., 2017}).

Looking back at the point P and point P’, they have the same plane position, but the difference is their positive height. Therefore, it can be seen that the elevation anomaly η is related to positive elevation H _g . Formula (3) shows that the change of positive elevation is equal to that of geodetic height. Therefore, the elevation anomaly η is related to geodetic height H.

Ternary polynomial regression function

If x, y and H are taken as three independent variables and elevation anomaly η is taken as dependent variable, then there are functions

1) Three parameters

2) Four parameters

3) Seven parameters

4) Eight parameters

5) Eleven parameters

By calculating the known data set (x, y, H, η), the corresponding polynomial parameters can be obtained and the functional relationship can be determined. Among them, three-parameter function requires at least three known point data, four-parameter function requires at least four known point data, seven-parameter function requires at least seven known point data, eight-parameter function requires at least eight known point data, and eleven-parameter function requires at least eleven known point data. If there is more redundant known data, it can be calculated according to the least square criterion (^{Liu et al., 2016}; ^{Wang et al., 2016}; ^{Wu & Miao, 2010}).

In the formula, a is a coefficient matrix composed of known data x, y, H and xy, P is a weight matrix, generally a unit matrix, L is an observation matrix composed of known data η, and X is the parameter matrix obtained.

Case analysis

In order to analyze the feasibility of the ternary polynomial regression fitting method proposed in this paper, the first 12 points are taken as the reference points and the last 7 points are taken as the points to be fitted by using the data provided in literature (^{Li, 2013}). Planar fitting and three-parameter ternary polynomial regression fitting are performed separately. The known data is shown in Table 1 and Table 2.

Table 1 Data of reference points (unit: m) 

Table 2 Data of points to be fitted (unit: m) 

After calculation, the parameters of plane fitting function are respectively α ₀ = 34.1044, α _x = -8.0703e-6, α ₂ = -3.1533e-5 and, and the parameters of three-parameter ternary polynomial fitting are respectively α ₀ = -3.6290e-7, a ₁ = -4.4750e-5, α ₂ = -2.3390e - 4 and. The fitting results are shown in Table 3.

Table 3 Plane fitting and three-parameter ternary polynomial regression fitting results. 

By using the above data, quadric surface fitting and seven-parameter ternary polynomial regression fitting are used respectively. The fitting results are shown in Table 4.

Table 4 Quadratic surface fitting and seven-parameter ternary polynomial regression fitting results. 

Using the above data, the ternary polynomial regression fitting of three parameters, four parameters, seven parameters, eight parameters and eleven parameters is carried out respectively. The fitting accuracy is shown in Table 5, and the trend of fitting accuracy is shown in Figure 2.

Figure 2 Ternary polynomial fitting accuracy trend 

Table 5 Ternary polynomial regression fitting accuracy 

From the results in Table 3, it can be seen that the smallest value of the three-parameter ternary polynomial regression fitting residual is 0.17 cm, the average residual is 0.89 cm, and the smallest value of the plane fitting residual is 0.23 cm, the average residual is 2.61 cm. It is concluded that the ternary polynomial regression fitting method has strong feasibility in the actual GPS elevation fitting.

From the results in Table 4, it can be seen that the fitting effect of seven-parameter ternary polynomial regression is better than that of three-parameter ternary polynomial regression. The average residual is reduced to 0.27 cm, and the overall residual is mostly within 1 cm. Similarly, the fitting effect of quadric surface is better than that of plane fitting. The minimum residual is reduced to 0.05 cm, and the average residual is reduced to 1 cm. It can be seen that the overall effect of seven-parameter ternary polynomial regression fitting is better than that of quadric surface fitting.

Conclusion

In this paper, in the process of GPS elevation fitting, the geodetic height is included in the calculation system, and the ternary polynomial regression method is used to solve the elevation anomaly. The method takes into account the influence of the geodetic height on the elevation anomaly and avoids the phenomenon of distortion caused by the omission of the geodetic height. The feasibility of the proposed method and its superiority in fitting the global elevation anomaly in a large range are verified by the experimental results. In addition, it can be seen from Table 5 that in the actual engineering application, the fitting effects of the seven-parameter and eight-parameter are the best. The internal fitting accuracy is less than 0.2 cm, and the external fitting accuracy is about 1.5 cm. This paper puts forward new ideas and directions for GPS elevation fitting. In order to solve practical problems better, it can be combined with other methods, and a series of problems related to accuracy need to be tested and discussed further.