Accurate computation of geoid-quasigeoid separation in mountainous region – A case study in Colorado with full extension to the experimental geoid region

1 Introduction

The geoid-quasigeoid separation (GQS) has traditionally used the Bouguer anomalies to approximate the difference between the mean gravity along the plumb line and the mean normal gravity (Heiskanen and Moritz 1967). This approximation is adequate in flat, low elevation areas, but it is not accurate in high and rugged mountains. Flury and Rummel (2009) added a so-called potential correction (Pot. Corr.) to the GQS and they showed that the correction can reach decimeters in deep valleys and high mountains. Sjӧberg (2010) added another second-order correction called the gravity gradient correction (G.G. Corr.). To distinguish the GQS with/without the higher order correction terms, the GQS using the Bouguer anomalies is called the simple GQS (sGQS), and the separation with the higher order correction terms is called the complete GQS (cGQS).

The rigorous formulation of the GQS is expanded into a Taylor series of the Bouguer disturbance in Appendix A. The formulas are equivalent to those in Sjӧberg and Bagherbandi (2017), even though they look slightly different. The differences in the equations appear in the use of the mean normal gravity and the normal gravity on the geoid, and the difference between the complete Bouguer disturbance and Bouguer anomaly.

Computation of the topographic effects on gravity and potential are split into near- and far-zone computations. It is shown in Appendix B that the first-order topographic effects of the far-zone are canceled out, with only the second and higher order effects left. This implies that the topographic effects of the far-zone contribution can be neglected.

Section 2 formulates the differences between the simple and complete GQS. Besides the potential and G.G. Corrs., the correction due to the difference between the gravity anomaly and disturbance is also added. Section 3 summarizes the information of the datasets used in the computations. Numerical computation methods are described in Section 4. The numerical results of the correction terms in Colorado are discussed in subsection 4.1. The results of the complete GQS for the whole experimental geoid 2020 (xGEOID20) region is described in subsection 4.2. Summary of the results and conclusion are given in Section 5. The computed complete GQS in this article was used in the Colorado experiment and xGEOID20 computations (Wang et al. 2021a, b).

2 Corrections to the sGQS

The sGQS has been computed by using the Bouguer anomaly Δg _B (Heiskanen and Moritz 1967, p. 327):

where N and ς are the geoid height and height anomaly, respectively, γ ¯ is the mean normal gravity, and H is the orthometric height.

For a better approximation, Flury and Rummel (2009) added a so-called potential difference correction [ V t ( Q ) − V t ( P ) ] γ ¯ − 1 , where V t ( P ) and V t ( Q ) are the gravitational potential of the topographic masses (topographic potential hereafter) at the point P on the Earth’s surface and Q on the geoid, respectively. For completeness, Sjӧberg (2010) added the G.G. Corr. − H 2 2 γ ¯ ∂ Δ g B ∂ H . These corrections are accurate to the second order of the series expansion of the GQS (Sjӧberg and Bagherbandi 2017). The complete GQS in terms of the Bouguer anomaly reads

The GQS in terms of the Bouguer disturbance is derived and given in Appendix A. To the power of H 2 , the complete GQS in terms of Bouguer gravity disturbance δ g B is

In this article, we use Equation (3) for the computation of the complete GQS. The difference between the simple (Equation (1)) and complete GQS (Equation (3)) is

The first, second, and third terms on the right side of Equation (4) are called the correction of gravity difference, potential difference, and gravity gradient, respectively. The physical meaning of the three correction terms is obvious: the gravity difference correction (G.D. Corr.) is caused by the differences between the gravity disturbance and gravity anomalies; the potential differences are caused by the differences between the topographic potential at the geoid and on the Earth’s surface; and the G.G. Corr. represents the contribution of the Bouguer disturbance gradient. All three corrections are evaluated in the Colorado region to show their magnitudes, or the errors if the sGQS is used in the geoid/quasigeoid conversion. Then, the complete GQS is computed for the xGEOID20.

3 Data used and computation methods

The following datasets were used:

1. The xG20DEM in 3″ grid spacing along latitude and longitude directions (Krcmaric 2022). This DEM combined TanDEM-X (Wessel et al. 2018), MERIT (Yamazaki et al. 2017), and the USGS 3D Elevation Program (3DEP) dataset (Lukas and Baez 2021). The RMS of the height differences between xG20DEM and 31,851 GPS on benchmarks was 2.1 m.

2. Terrestrial gravity data consisting of 1,633,376 point-gravity data provided by the National Geospatial-Intelligence Agency and 135,290 NGS point gravity data over the ocean areas, and the DTU15 (Andersen et al. 2016) altimetric gravity anomalies in 1 ' × 1 ' grid spacing were used.

3. The xGEOID19B (Li et al. 2019) was used for computing the G.D. Corr. for the whole xGEOID20 region. The maximum geoid differences between the GEOID18 (Ahlgren et al. 2020) and xGEOID19B over the contiguous United States was about a meter, which caused a gravity error of 0.3 mGal in the gravity anomaly differences, and hence ignored.

The simple and complete GQS were computed using the abovementioned data-sets. Note that the gravity dataset used in the Colorado experiment was a subset of the terrestrial gravity (Wang et al. 2021a). However, the DEM used for the experiment was the SRTM v4.1 (Jarvis et al. 2008) and the updated xG20DEM was used in this article.

4 Numerical computation and results

Computations of the topographic effects on gravity and potentials were split into near-zone and far-zone contributions. It is shown in Appendix B that the first-order topographic effect on the GQS is canceled out in the far-zone contribution, only second and higher order terms remain. Therefore, the far-zone contribution to the GQS was neglected in this article.

The terrain correction was computed at every location of the terrestrial gravity observations using the xG20DEM elevation model in 3″ × 3″ grid size, and a one-degree spherical cap (near-zone) was used for the integration area. Then, the point Bouguer anomaly was computed and gridded into 1′ × 1′ grid using the program GEOGRID (Forsberg and Tscherning 2008). The Bouguer anomaly is then analytically (harmonically) downward continued onto the geoid using the integral equations (8–85) in Heiskanen and Moritz (1967). The G.G. Corr. is then computed from the harmonic downward continued Bouguer anomaly using equations (8–87) in (ibid.) to the second approximation.

The difference between the Bouguer gravity disturbance and Bouguer gravity anomaly can be computed as

where Δg and δ g are the gravity anomaly and gravity disturbance, respectively.

4.1 Results for Colorado

Because the correction terms have meaningful contribution only in high and rugged mountains, the Colorado region (35° ≤ ϕ ≤ 40°, −110° ≤ λ ≤ 102°) was selected as a test area. There is another reason for this selection: the 1 cm-geoid experiment was conducted in this area and the computed GQS was used to convert the geoid models into the quasigeoid models for comparisons and evaluations (Wang et al. 2021a). Three correction terms are shown in this section. To show the relationship between the corrections and the topography, the topographic height is plotted in Figure 1.

Figure 1

Topographic height of xG20DEM in the Colorado area, spatial resolution 3″. The elevation has a mean value of 2,014 m and a STD DEV value of 614 m. The minimum height is 929 m and the highest peak is 4,381 m. The red line represents the GSVS17 (van Westrum et al. 2021) profile along the highway US 160 from Durango to Walsenburg (west to east).

The xGEOID19B ranges from −26.98 to −12.53 m of geoid height in this region. Thus, the gravity differences computed from equation (5) are negative with a mean value of −19.43 m resulting in a negative mean value of −6.0 mGal in the gravity differences. The contribution of the gravity differences to the GQS is plotted in Figure 2. For comparison, the same color scale is used in Figures 2–5.

Figure 2

Gravity difference correction to the GQS. Units in meters. Mean value = −0.011, Min = −0.020, Max = −0.007, Std Dev = 0.002, and RMS = 0.012.

Figure 3

Potential correction. Unit in meters. Mean value = −0.004, Min = −0.190, Max = 0.117, STD DEV = 0.018, and RMS = 0.019.

Figure 4

Contribution of the gravity gradients correction to the separation. Unit in meters. Mean value = 0.000, Min = −0.025, Max = 0.020, Std Dev = 0.003, and RMS = 0.003.

Figure 5

Difference between the complete and sGQS. Units in meters. Mean value = −0.08, Min = −0.135, Max = 0.180, Std Dev = 0.017, and RMS = 0.019.

Figure 2 shows that the gravity correction is small – merely 0.002 m in terms of Std Dev values. However, there is a −0.011 m bias and its extreme value reaches −0.020 m, too. As this correction is easy to compute and it reaches cm-level, it can be easily included in the GQS computations.

The topographic potential is computed at the Earth’s surface and on the geoid using a 1′ × 1′ grid with a one-degree spherical cap as the integration area (near-zone). As shown in Appendix B, the far-zone contribution is in second and higher order terms and can be ignored. Figure 3 illustrates the Pot. Corr. to the GQS.

Figure 3 shows that the Pot. Corr. has decimeter contribution to the GQS. The correction ranges from −0.190 to 0.117 m, with a standard deviation of 0.018 m. This correction is the largest of all three corrections. Another observation from Figure 3 is that the correction shows strong correlation with the terrain variation.

The G.G. Corr. is a second-order term and it is plotted in Figure 4.

Similar to the G.D. Corr., the G.G. Corr. is small in terms of Std Dev and RMS values, and it ranges from −0.025 m to 0.020 m. In other words, the correction of gravity gradients is statistically small, but it can have a few cm contributions pointwise.

The combination of all three corrections, or the difference between the complete and sGQS, is plotted in Figure 5.

From all the above results, we have the following conclusions:

1. The largest correction comes from the potential difference term in Equation (4). Its extreme values reach −0.19 m, with an RMS value of 0.019 m.

2. Contributions of the other two corrections are small. The G.G. Corr. is about 0.003 m in terms of RMS value, but it can reach 0.03 m pointwise in the extreme. The gravity correction (free-air vs gravity disturbance) has a −0.011 m bias and its extreme values reach −0.02 m because of the high elevation of the region.

In summary, the contribution of the correction terms to the GQS reaches 0.180 m in the Colorado area of steep topography where the average height of topography is 2,017 m and the highest peak reaches 4,385 m, and the RMS of the corrections is 0.019 m.

In the end, it is necessary to point out that the Bouguer anomaly is the main term in the GQS: it is an order larger than the correction terms. The simple separation (using the Bouguer anomalies only) is plotted in Figure 6.

Figure 6

The sGQS. Units in meters. Mean value = −0.458, Min = −1.376, Max = −0.126, Std Dev = 0.230, and RMS = 0.513.

As a comparison, the cGQS in Colorado is plotted in Figure 7.

Figure 7

The complete GQS. Units in meters. Mean value = −0.454, Min = −1.257, Max = −0.127, Std Dev = 0.223, and RMS = 0.506.

Comparing Figure 7 to Figure 6, the complete separation is slightly smoother than the simple separation, evidenced by smaller extreme values (−1.257 vs −1.376 m) and smaller Std Dev (0.223 vs 0.230 m).

4.2 Results for the xGEOID20 region

After the regional runs for Colorado, the GQS was computed for the xGEOID20 region. Over ocean areas, the elevation is zero, therefore making the separation zero. Thus, the statistical analysis was on land only. To have a full picture, the statistics of the sGQS, cGQS and the correction terms are listed in (Table 1).

Table 1

Statistics of the simple separation (sGQS) and complete separation (cGQS), potential correction (Pot. Corr.), gravity gradient correction (G.G. Corr.), and gravity difference correction (G.D. Corr.), land only. Number of samples: 15,252,799. Units are in meters

	sGQS	cGQS	Pot. Corr.	G.G. Corr.	G.D. Corr.
Mean value	−0.057	−0.055	0.001	−0.000	0.002
Std Dev	0.119	0.119	0.006	0.002	0.011
RMS	0.132	0.131	0.006	0.002	0.011
Min	−1.348	−1.275	−0.144	−0.090	−0.031
Max	1.082	1.296	0.590	0.102	0.065

Statistically, the sGQS and cGQS are very close. The corrections are around a cm or less in terms of the RMS, but the corrections can reach 0.590 and 0.102 m pointwise for the potential and the G.G. Corrs., respectively. Large corrections happen in high mountains (Table 2) with rugged terrain. Another observation from Table 1 is that the mean value of the cGQS is −0.055 m, showing that the mean value of the height anomalies is 5.5 cm higher than the mean value of the geoid undulation over land areas. Note that the height anomalies and geoid heights are the same over ocean areas where the GQS is zero.

Table 2

RMS values of the GQS as a function of height. Units in meters

Elevation (m)	No. of pts	sGQS	cGQS	Pot. Corr.	G.G. Corr.	G.D. Corr.
0–500	8,778,198	0.013	0.015	0.001	0.000	0.002
500–1,000	2,544,283	0.055	0.058	0.003	0.001	0.005
1,000–1,500	1,376,066	0.143	0.147	0.007	0.002	0.009
1,500–2,000	1,051,370	0.236	0.239	0.010	0.002	0.015
2,000–2,500	807,425	0.318	0.315	0.012	0.003	0.023
2,500–3,000	553,678	0.293	0.285	0.014	0.005	0.034
3,000–3,500	132,939	0.513	0.490	0.025	0.010	0.040
3,500–4,000	7,792	0.754	0.717	0.088	0.008	0.023
4,000–4,500	806	0.451	0.387	0.149	0.010	0.023
4,500–5,000	187	0.481	0.431	0.254	0.013	0.027
5,000–5,500	48	0.423	0.293	0.374	0.016	0.032
5,500–6,000	7	0.433	0.136	0.522	0.016	0.033
0–6000	15,252,799	0.132	13.1	0.006	0.002	0.011

It has been assumed that the magnitude of GQS would increase with height, and the largest correction would be found around the highest mountains. Statistical analysis of the computed GQS supports this assumption partially, but with some modification. Table 2 lists the RMS values of the separation with respect to topographic height.

Table 2 shows that the RMS of GQS (simple and complete) increase with height until 4,000 m, then starts decreasing. Unexpectedly small RMS of cGQS happens in seven cells with the height between 5,500 and 6,000 m, the highest peaks of the xGEOID20 region. Another observation is that the simple and complete GQS terms have similar RMS values until 3,000 m, then start to differ. The largest difference (0.433 vs 0.136 m) happens for seven cells with heights between 5,500 and 6,000 m. But the reader is cautioned for the interpretation of the statistics at 5,000 m and higher elevations because of its small number of samples, which are not statistically significant.

To inspect the unexpected small GQS in detail, the Bouguer anomaly (BA), gravity gradient of the Bouguer disturbance G1, the terrain correction (TC) and the GQS terms for seven cells are listed in Table 3.

Table 3

Gravities and GQS in cells higher than 5,500 m

Mount Logan
Latitude (degree)	Longitude (degree)	Orth H (m)	BA (mGal)	G1 (mGal)	TC (mGal)	sGQS (m)	cGQS (m)	Pot. Corr. (m)	G.G. Corr. (m)	G.D. Corr. (m)
60.56667	−140.4000	5,589	−79.9	8.1	146.1	−0.455	0.150	0.548	0.023	0.034
60.58333	−140.4500	5,531	−76.1	9.5	106.6	−0.428	0.105	0.472	0.027	0.034
60.58333	−140.4333	5,624	−78.6	8.7	112.0	−0.450	0.118	0.508	0.025	0.035

Mountain Denali
Latitude (degree)	Longitude (degree)	Orth H (m)	BA (mGal)	G1 (mGal)	TC (mGal)	sGQS (m)	cGQS (m)	Pot. Corr. (m)	G.G. Corr. (m)	G.D. Corr. (m)
63.06667	−151.0167	5,715	−73.9	0.7	157.4	−0.430	0.173	0.570	0.002	0.031
63.06667	−151.0000	5,768	−74.0	0.7	158.8	−0.435	0.189	0.590	0.002	0.032
63.08333	−151.0333	5,549	−73.6	0.5	129.6	−0.416	0.093	0.477	0.001	0.030
63.08333	−151.0167	5,589	−73.7	0.5	123.6	−0.420	0.088	0.475	0.001	0.031

The steep slopes around the peaks cause large terrain corrections and large Pot. Corrs. that compensate for the contribution of the Bouguer disturbances, resulting in small GQS. The Pot. Corr. is around a half meter for those cells, which has the same magnitude as that of the sGQS. The G.G. Corr. is in mm level for the peaks of Denali, but over 2 cm at the peaks of Mount Logan. The G.D. Corrs. at the peaks are around 3 cm.

The GQS is a combination effect of gravity anomalies, the terrain roughness, and the height of the topography. Thus, the extreme values of GQS happen not necessarily at the highest mountains. The information for the extreme values is listed in Table 4.

Table 4

Information for the maximum and minimum values of the GQS separation

Latitude (degree)	Longitude (degree)	Orth H (m)	BA (mGal)	G1 (mGal)	TC (mGal)	sGQS (m)	cGQS (m)	Pot. Corr. (m)	G.G. Corr. (m)	G.D. Corr. (m)
19.81667	−155.4667	3,992	265.2	9.7	45.6	1.082	1.296	0.161	0.020	0.033
37.85000	−107.4500	3,756	−331.3	−10.0	2.7	−1.270	−1.275	0.032	−0.019	−0.018

The maximum separation happens at Mauna Kea, Hawaii. The cell has an elevation of 3,992 m (xG20DEM), and the Bouguer anomaly is positive 265 mGals. The Pot. Corr. to the GQS reaches 16 cm. The minimum separation is at Uncompahgre Peak, Colorado with the height of 3,756 m (xG20DEM) where the terrain seems moderate, indicated by a small terrain correction of 2.7 mGal (Table 4).

Figure 8

The cGQS in the experimental geoid region. Units in meters. Decimeter to meter separations happen in the rocky mountains. The large corrections in Greenland are due to sparse surface gravity observations, resulting in large errors in the Bouguer anomalies. Improvement will be made in the future by using the ARCGP gravity and the ice thickness data in the area.

5 Summary and conclusion

Three correction terms were added to the simple GQS based on the Bouguer gravity anomalies. The correction terms are first tested in Colorado. The results show that the corrections are small in terms of RMS value (0.019 m), but reach a few cm to decimeters (0.19 m) pointwise. The corrections are significant for the cm geoid computations and must be included in computations.

The cGQS was then computed for xGEOID20. Statistically, more than half of the 1′ × 1′ cells are below 500 m and the RMS of the GQS is merely 0.015 m. The results also support the general assumption that the magnitude of the GQS increases with respect to the topographic height, but only to elevations of 4,000 m where the complete GQS reaches its maximum of 0.717 m in terms of RMS. The RMS values then decrease with respect to the elevations. The RMS values are reduced to 0.293 m at heights of 5,500 m and 0.136 m at heights of 6,000 m. This surprising occurrence mostly is attributed to the steep slopes around the peaks that contribute large Pot. Corrs. (about half meter), and it compensates for the contribution of the Bouguer anomaly and results in a smaller GQS.

The magnitude of GQS is a combination effect of gravity anomalies, terrain roughness, and height of the topography (Figure 8). The maximum GQS 1.296 m happens not at the highest peaks or the deepest valleys, but at Mauna Kea, Hawaii where the elevation is 3,993 m (xG20DEM), and the Bouguer anomaly is positive 265 mGals. The minimum GQS of −1.275 m is located at the Uncompahgre Peak, Colorado, which has a height of 3,756 m.

The far-zone topographic effects are in the second and higher order terms, which are neglected in this computation. The higher order of vertical derivatives of Bouguer disturbance is also omitted. Even if they are considered as negligible for most cases, the effects will be addressed in the next experimental geoid computations to make sure that they are below 1-cm pointwise.