Evaluating the Impact of the Recent Combined and Satellite-Only Global Geopotential Model on the Gravimetric Geoid Model

Geoid represents Earth’s surface, ocean, and gravitational field, which influence the elevations, shape, and mass distribution of the geopotential surface, a hypothetical surface that is perpendicular to the direction of gravity at every point. This geopotential surface serves as a reference for measuring elevations and is used as a fundamental reference surface for geodetic and surveying purposes. In this study, the Least Squares Modification of Stokes Formula (LSMS) with Additive Corrections (AC), also known as the KTH method, is used to generate a new gravimetric geoid model for Peninsular Malaysia. The KTH method was developed at the Royal Institute of Technology (KTH) in Stockholm-Sweden. The dataset used is the most recent global digital elevation model, Shuttle Radar Topography Mission (SRTM) 1 Arc-Second Global, generated by the National Aeronautics and Space Administration (NASA) and the National Imagery and Mapping Agency (NIMA). In addition to this elevation data, the dataset includes the Global Geopotential Model (GGM), which is composed of the XGM2016, XGM2019e, Tongji_GGMG2021S, and Tongji-Grace02k models. Furthermore, it incorporates sets of regional gravity data, including terrestrial gravity, airborne gravity, and marine gravity anomalies, all of which are derived from the Technical University of Denmark (DTU 21). The actual 45 Global Navigation Satellite System (GNSS)-levelling points data have been compared to the gravimetric geoid model developed in this study and the geoid acquired from Department of Survey and Mapping Malaysia (DSMM). According to the statistical results, NXGM2019e provides better accuracy, with the Root Mean Square Error (RMSE) geoid model errors of ±0.033 m, compared to the deviations in free-air anomalies, XGM2019e, which has the minimum RMSE of 10.291 mGal. Meanwhile, Tongji-GMMG2021S has the maximum RMSE of 14.792 mGal. The geoid is derived from the XGM2019e model and has maximum and minimum values of 0.032 m and 0.147 m, respectively, with mean residuals of 0.089 m. In conclusion, the XGM2019e has the potential to determine a precise local geoid model for Peninsular Malaysia


Introduction
The determination of a precise gravimetric geoid model is of utmost importance in geodesy and geophysics.Geoid serves as a fundamental benchmark in determining elevations in relation to mean sea level and holds significance in a range of practical contexts, including surveying, mapping, navigation, and geophysical studies [1].It provides a consistent reference surface for surveying and mapping, enabling elevation measurement, vertical datum establishment, and standardisation of height measurements across different locations.Geoid also serves as the basis of important decisions in engineering, construction, and scientific studies.The precision of a gravimetric geoid model is contingent on various factors, with the calibre and extent of the gravity field data employed in its calculation being particularly important.To overcome the limitations of local gravity datasets, geoid calculation involves long-wavelength data from Global Geopotential Model (GGMs) with local gravity data obtained from terrestrial measurements, following Stokes' principles.Nova k et al. [2] demonstrated that parameters that characterise the external gravitational field can be derived from the inversion of satellite data.
Satellite gravity missions, such as Gravity Recovery and Climate Experiment (GRACE) [3], CHAllenging Minisatellite Payload (CHAMP) [4], and Gravity Field and steady-state Ocean Circulation Explore (GOCE) [5], have provided significant contributions to the improvement of global gravity field models.However, the effectiveness of these models in specific regions depends on their inherent quality and parameter settings [6].Several studies by Doganalp [7], Nicacio et al. [8], Goyal et al. [9], Lee and Kwon [10], and Nyoka et al. [11], have conducted evaluations in their respective fields to identify the most suitable model.Recent research by Wu et al. [12] compared two gravity field models, XGM2016 and EGM2008, and found that the selection of the gravity field model significantly affects the accuracy of the geoid model within the study area.This finding is consistent with the results of a study conducted by Isik et al. [13], which involved the construction of a gravimetric geoid model for the Turkey region using several satellite-only models.Hence, the meticulous selection of a suitable gravity field model dataset is imperative to attain a precise geoid model.
The main objective of this research is to examine the influence of the XGM and Tongji models on gravimetric geoid modelling in Peninsular Malaysia.To achieve this objective, four geoid models have been computed using four different GGMs: XGM2019e, XGM2016, TONGJIGGMG2021S, and Tongji-Grace02k.The geoid models presented in this study have been derived using the methodology established by Sjoberg in 1984Sjoberg in , 1991Sjoberg in , 2003aSjoberg in , and 2003b, as employed by the Royal Institute of Technology (KTH).The KTH method differs from other approaches which directly incorporate gravity anomalies without reduction [14].Instead, it modifies an initial geoid approximation to minimise modelling errors [15].This method determines least-square parameters to minimise overall errors, thereby mitigating the impact of geoid modelling inaccuracies [16].
The KTH method has been utilised in the computation of gravimetric geoid models in several regions and countries, including Sudan [17], Malaysia [18], Poland [19], Uganda [20], Bosnia and Herzegovina [21], and Latvia [22].These countries are among the prominent users of the KTH method in this field.The first section of this paper explains the introduction and concept of the research.The second section of this paper describes the input data, while the third section discusses the mathematical aspects of the KTH method and the additive corrections.The fourth and fifth sections present the novel gravimetric geoid model (KTH-SDG08) and its validation, respectively.Finally, concluding remarks are provided in the sixth section.

Study Area
The geographical region for geoid computation in this research is Peninsular Malaysia, with the optimal cap size at 0≤lat≤8.75 and 97≤long≤106.75as illustrated in figure 1.

Terrestrial gravity data
The gravimetric geoid model for Peninsular Malaysia has been developed by utilising terrestrial gravity data obtained from the Department of Survey and Mapping Malaysia (DSMM).The dataset, illustrated in Figure 2a, comprises 8,500 gravity points covering the entire study area.In cases where terrestrial gravity data are sourced from multiple providers, it often contains redundant and erroneous points that necessitate identification and removal.To address these concerns within the gravity database, this study employs a two-step data cleansing procedure.Initially, duplicate sites are identified by comparing their longitude, latitude, and gravity values.The second phase involves the implementation of a cross-validation approach, which is extensively described by Yilmaz et al. [23].

Airborne gravity data
The geoid computation in this study relies on a secondary gravity dataset obtained from DSMM, which consists of airborne gravity data.This dataset is particularly valuable in regions with limited terrestrial gravity data, such as mountainous and forested areas.A total of 24,855 airborne gravity measurements have been utilised for the geoid computation, and the spatial distribution of this dataset is illustrated in Figure 2b.It is important to note that the preprocessing of the raw airborne gravity data has been conducted by DSMM prior to this study, therefore it is not within the scope of this research.

Marine gravity data
The global marine gravity model utilised in this study was sourced from the latest version, DTU21, which was obtained from the DTU National Space Institute in Denmark.To overcome the significant difficulties associated with gravity anomalies derived from satellite altimetry in nearcoastal marine regions, a specific exclusion zone has been implemented.This exclusion zone extended 20 km from the coastline and 15 km from any airborne gravity data, ensuring that the marine gravity anomalies within the model are not influenced by these areas.The geographical distribution of the global marine gravity model, DTU21, is depicted in Figure 2c.

Digital Elevation Model (DEM)
To accurately correct terrain and Bouguer effects in the development of a gravimetric geoid, the use of the Digital Elevation Model (DEM) is crucial.This study utilises a high-resolution Digital Elevation Model (DEM) to ensure the utmost precision in geoid computation.The study utilises DEM data acquired from the Shuttle Radar Topography Mission (SRTM), which has a spatial resolution of 1 arc-second, corresponding to a ground resolution of approximately 30 meters.This data source, depicted in Figure 2d, has been derived from a collaborative mission conducted by the National Aeronautics and Space Administration (NASA).

GNSS/Levelling data
The dataset employed in this research comprises GNSS levelling data, encompassing the ellipsoidal height (h) derived from GNSS observations and the corresponding benchmark levelling height (Hmsl).This dataset is crucial to ensure the accuracy of geoid models and to evaluate the impact of different density models used in geoid modelling.Two distinct sets of GNSS levelling data have been employed in this study.The initial dataset, supplied by DSMM, comprises of 45 sets of coordinates.These coordinates play a vital role in assessing the precision of the geoid model throughout the entire area of Peninsular Malaysia.The main objective of these measurements is to investigate the relationship between mountainous terrain and the performance of the geoid model.The spatial distribution of the initial set of coordinates is illustrated in Figure 3.

Global Geopotential Model
Gravity field models can be categorised into three primaries: satellite-only, combined, and tailored gravity field models.The International Centre for Global Earth Models (ICGEM) website provides public access to approximately 178 GGMs (http://icgem.gfz-potsdam.de/home).One of the latest additions is the Tongji-GMMG2021S, released in 2022 as a satellite-only GGM.Unlike its predecessor, Tongji-Grace02k, which is solely based on GRACE satellite data, Tongji-GMMG2021S combines GRACE and GOCE data.This combination involves reprocessing GOCE Level 1b gravity gradient observations and incorporating with the 180 degrees of freedom (d/o) model Tongji-Grace02s.These models represent mathematical representations of the Earth's gravitational field, utilising data collected from GRACE and GRACE Follow-On (GRACE-FO) satellite missions.
XGM2019e is a combined global gravity field model that incorporates spheroidal harmonics up to degree and order (d/o) 5399 and has a spatial resolution of approximately 2 minutes (roughly 4 kilometres).Released in 2019, this GGM integrates data from three primary sources: the satellite model GOCO06s, which contributes data up to d/o 300 in the longer wavelength range, ground gravity measurements, and ocean gravity information.The ground gravity data includes a 15-minute ground gravity anomaly dataset provided by the US National Geospatial- Intelligence Agency (NGA), supplemented by topographically derived gravity data over land (EARTH 2014).For the ocean region, gravity anomalies from DTU13 with a 1-minute resolution have been utilised.More detailed information about the XGM2019e model can be found in Zingerle et al. [16].Notably, XGM2019e represents a significant improvement over its predecessor, EGM2008, in terms of its combined models with a maximum degree/order exceeding 2000, as documented by the ICGEM.Unlike previous higher-resolution combined GGMs that relied on EGM2008 as a data source, XGM2019e operates independently.
Several studies have investigated the impact of the combined global gravity field model XGM2019e on gravimetric geoid modelling.In a recent study conducted in Nigeria, the XGM2019e model was employed to calculate gravity field parameters.The purpose of the study was to assess the effectiveness of different combined Global Geopotential Models (GGMs) for geological mapping [24].In Kuwait, a geoid model with high resolution was calculated using XGM2019e as a component of the combined global geopotential model [25].Furthermore, Vergos et al. [26] employed XGM2019e as a benchmark for simulating the long-wavelength segment of the gravity field spectrum in a localised gravimetric geoid model.The XGM2016 global geopotential model offers comprehensive insights into the gravitational field of the Earth, employing satellite data obtained from missions including the Gravity Recovery and Climate Experiment (GRACE) and the Gravity field and steady-state Ocean Circulation Explorer (GOCE).These satellites measure gravitational variations caused by Earth's mass distribution, enabling the development of models like XGM2016.XGM2016 demonstrates higher accuracy compared to EGM2008, both over land and ocean.
Detailed discussions comparing XGM2019e_2159 with XGM2016 can be found in various sources.Zingerle et al. [16] noted that XGM2019e_2159 performs exceptionally well over oceans, showing comparable or even superior performance compared to earlier models.Furthermore, Wu et al. [27] highlighted that alternative Global Geopotential Models (GGMs) exhibit reduced accuracy compared to XGM2019e_2159, with coastal zones being significantly affected.Table 1 presents data regarding the content of satellite data and the highest degree of order for the four models.Tongji-GGMG2021s 2022 300 S (Goce), S (Grace) [30] Tongji-Grace02k 2018 180 S (Grace) [30] a where A is for altimetry, S is for satellite (e.g., GRACE, GOCE, LAGEOS), G is for ground data (e.g.terrestrial.shipborne, and airborne measurements) and T is for topography.

Accuracy assessment of GGM
Prior to gravimetric geoid computation, a thorough evaluation of all tested Global Geopotential Models (GGMs) has been conducted.This assessment utilises the data collected from 45 GNSS levelling points and 8000 terrestrial gravity points.The assessment of GGMs involves the computation of geoid height and free air anomalies obtained from these models.The values obtained have been subsequently interpolated onto the corresponding positions of GNSS levelling and terrestrial gravity points using the computation service offered on the ICGEM website.Throughout this computation, a Tide Free system has been consistently applied for the Tide system, and the GRS80 reference system has been used.This approach ensures the consistency of the comparison and reduces potential biases arising from differences in reference ellipsoids, tide conventions, and geoid references.
To evaluate the performance of each GGM, a fundamental statistical analysis has been conducted.This analysis involves calculating the Mean Absolute Error (MAE) and Root-Mean-Square Error (RMSE) by comparing the geoid height obtained from the GGMs (denoted as ), along with geometric geoid height derived from GNSS levelling (denoted as  ).The calculations have been carried out as follows: where  represents the height or elevation of a point above the reference ellipsoid. represents the height or elevation of the same point above the mean sea level.Δ represents the difference in geoid undulation between two points.
The MAE and RMSE can be computed using the equations ( 3) and (4), respectively.
The assessment of GNSS levelling and geoid heights obtained from GGMs can be influenced by several sources of error.The errors discussed in this study include inconsistencies in data, as investigated by Yilmaz et al. [23] and Goyal et al. [9], as well as inaccuracies that are inherent in the GGMs, such as commission and omission errors.Furthermore, the GNSS and levelling procedures may introduce biases and errors [31].To minimise the effects of systematic biases, the discrepancies have been corrected using a 4-parameter model that is founded on: = [1 cos   cos   cos   sin   sin   ]  (6) In this framework, the vector  represents a set of unknown parameters, whereas  stands for known coefficients.Additionally, the free air anomalies, denoted as Δ, are determined from terrestrial gravity data before conducting a comparison with the GGMs and computing the gravimetric geoid model.This computation has been executed utilising equations ( 8) and ( 9), which are formulated as follows: In this context,  represents the gravity value at the topographic surface,  denotes the correction for the free air (approximately 0.3086 times ), and  stands for the normal gravity on the GRS80 ellipsoidal surface.The  is derived using the Somigliana Formula, and it is expressed as follows: In this equation,  represents gravity at the equator,  is gravity at the pole,  signifies the semimajor axis of the ellipsoid,  represents the semi-minor axis of the ellipsoid,  denotes the first eccentricity, and  stands for the geodetic latitude.The statistical analysis of the difference in gravity anomalies, Δ, between the terrestrial gravity-derived free-air anomalies, Δ, and GGMs is computed utilising Equation (8).Further details on the methodology used to evaluate GGMs are provided in Figure 4.

Geoid Computation utilising the KTH method
Before developing geoid using the KTH method, it is crucial to integrate all available gravity datasets, including terrestrial, marine, and airborne data, onto a uniform grid with a specified resolution.In this study, a grid spacing of 1'x1' has been employed to grid the gravity data, resulting in a gravimetric geoid with the same grid spacing.Several approaches for combining and gridding gravity data have been discussed in previous studies by Featherstone and Kirby [32], Goos et al. [33], McCubbine et al. [34], and Pa'suya et al. [15].However, for this study, the approach proposed by McCubbine et al. [34] and Pa'suya et al. [15] has been implemented.Detailed information regarding this strategy can be found in both studies.In the initial phase of geoid calculation using the KTH method, the gridded surface gravity anomalies, and the Global Geopotential Model (GGM) have been utilised to estimate the approximate geoid undulation (̃ ) through Stokes integration, as expressed below [35]: In the given equation, 0 represents the spherical cap,  stands for the mean Earth radius,  signifies the mean normal gravity,  represents the unreduced surface gravity anomaly, () denotes the modified Stokes's function,  represents the maximum degree of the Global Geopotential Model (GGM), and ∆  represents the Laplace surface harmonic of the gravity anomaly determined by the GGM of degree .Within the framework of the KTH method, four additive corrections are calculated, which include the combined topographic correction ̃  , the downward continuation correction ̃, the total atmospheric correction ̃, and the total ellipsoidal correction ̃.These corrections are subsequently combined with the approximate gravimetric geoid in the following manner:  =  ̃+   +   +   +   (12) The initial additive correction, referred to as the combined topographic correction, is calculated using the formula presented in reference [36].
where G is the Newtonian gravitational constant, ρ is the Earth's crust density, R is the Earth's radius, and H is the elevation of the topography at the computation point P. The second correction is DWC, which is calculated using Equations ( 14) -( 18): =   (1) In this context, the symbol  denotes the spherical radius of point P, where = R + , with  representing the orthometric height of point P. Additionally, Q is utilised as the moving integration point.The third correction, referred to as the simple ellipsoidal correction by Ellmann and Sjo berg [37], is computed as follows: The ultimate correction, as proposed by Sjo berg [36], involves a combined atmospheric correction: This equation incorporates , which represents the zero-degree term of atmospheric density at sea level, and , which signifies the Laplace surface harmonic of degree  pertaining to the topographic height.A comprehensive representation consisting of combined gravity data employed in this study is provided in Figure 5.

Evaluation of the GGMs performance
As outlined in Section 3.1, the performance of the four Global Geopotential Models (GGMs) has been assessed using terrestrial gravity-derived free-air anomalies and GNSS levelling-based geometrical geoid data.The findings of these evaluations are summarized in Tables 2 and Table 3, which present statistical analyses for terrestrial gravity anomalies and GNSS-levelling geometrical geoid heights, respectively.Generally, the comparison with terrestrial gravity anomalies indicates that the satellite-only GGMs, Tongji-Grace02k and Tongji-GMMG2021S, exhibit higher precision compared to the combined GGMs.This observation is in line with expectations, as combined GGMs are susceptible to contamination from errors with long and medium wavelengths in terrestrial gravity data.Among the satellite-only GGMs, Tongji-Grace02k fit with terrestrial gravity anomalies, with a mean error (ME) of 19.491 mGal and an RMSE of 14.137 mGal.This performance surpasses Tongji-GMMG2021S, which has ME and RMSE values of 19.440 mGal and 14.792 mGal, respectively.This disparity is somehow surprising, given that the Tongji-GMMG2021S model represents the most recent iteration of the global geopotential model available on the ICGEM website.This model is constructed by integrating data from GRACE and GOCE.In the case of the combined GGMs, both XGM models demonstrate nearly identical levels of accuracy.Although XGM2019e exhibits slightly better accuracy than XGM2016, with ME and RMSE values of 19.475 mGal and 10.291 mGal, respectively.This difference is not significant when compared to XGM2016, which presents accuracy figures of 19.462 mGal for ME and 10.354 mGal for RMSE, respectively.
Table 3 presents statistical parameters that describe the discrepancies between geoid undulations obtained from GNSS/levelling data and those derived from four different Global Geopotential Models (GGMs): XGM2016, XGM2019e, Tongji-GMMG2021S, and Tongji-Grace02k.In general, the accuracy of the four Geoid Height Models (GDEMs) is relatively similar.However, the combined GGMs tend to show better agreement with GNSS levelling data in comparison to the satellite-only models.Among the combined GGMs, XGM2019e demonstrates slightly superior accuracy relative to XGM2016, with ME and RMSE of 0.094 m and 0.306 m, respectively.On the other hand, among the satellite-only GGMs, the Tongji-Grace02k model once again indicates a better fit with GNSS levelling data compared to the latest model, Tongji-GMMG2021S, with ME and RMSE values of 0.089 m and 0.298 m, respectively.These comparisons reveal that distinguishing the best fit GGM is challenging due to the subtle distinctions in accuracy among them.

Evaluation of the gravimetric geoid derived based on difference GGM
The gravimetric geoid models computed using four different global geopotential model have been assessed by comparing to the GNSS-levelling points used in the previous evaluation of GGMs.To facilitate the comparison, the gridded gravimetric geoid has been interpolated to the positions of the GNSS-levelling points using bi-cubic interpolation.The geoid heights are then compared to the geometric geoid heights obtained from the GNSS-levelling data.As previously stated in Section 3.1, the discrepancies observed between the gravimetric geoid height and the geometric geoid height are affected by a range of systematic biases.These differences have been subsequently fitted using the 4parameter model, and the standard deviation of the residual has been computed.Figure 4 illustrates the gravimetric geoid models obtained from the four different GGMs, while Table 4 provides the statistical analysis comparing each geoid model with the GNSS-levelling data.The results indicate that the gravimetric geoid derived from XGM2019e exhibits the closest alignment with the GNSS-levelling data, with an RMSE of 0.033 meters when compared to XGM2016.The precision of the gravimetric geoid obtained from XGM2016 exhibits only marginal disparities compared to XGM2019e, with an approximate accuracy of 0.035 meters.Interestingly, the precision of the gravimetric geoid derived from satellite-only GGMs closely matches that of the combined GGMs, as shown in Table 4.Among the satellite-only GGMs, the gravimetric geoid model Tongji-Grace02k demonstrates higher accuracy compared to Tongji-GMMG2021S, with an RMSE of 0.035 meters.The gravimetric geoid derived from Tongji-GMMG2021S exhibits an RMSE of 0.041 meters.
In terms of ME, the gravimetric geoid derived from Tongji-Grace02k exhibits the lowest ME, followed by Tongji-GMMG2021S, XGM2016, and finally XGM2019e, with corresponding ME values of 0.071 meters, 0.089 meters, 0.094 meters, and 0.095 meters.This indicates that the selection of the Global Geopotential Model (GGM) significantly influences the accuracy of the gravimetric geoid model.Consequently, the selection of the most suitable GGM prior to geoid modelling is of paramount importance.

Conclusion
This research aims to evaluate the impact of a combined global geopotential model on the development of a gravimetric geoid model for Peninsular Malaysia.The study utilises the Least Squares Modification of the Stokes Formula (LSMS) with Additive Corrections (AC), also known as the KTH method, to generate the geoid model.Various data sources, including terrestrial gravity data, airborne gravity data, marine gravity data, DEM, and GNSS levelling data have been employed in the research.The primary objective is to assess the accuracy, spatial resolution, and consistency of the newly generated geoid model.The KTH method involves several sequential steps, including the application of Stokes' formula, interpolation of gravity data, incorporation of error degree variance, computation of modification parameters, and determination of the geoid model.The results and analysis section presents the comparisons of the geoid models derived from different GGMs, namely XGM2016, XGM2019e, Tongji-GMMG2021S, and Tongji-Grace02k.Statistical measures, including minimum, maximum, mean, and RMSE have been conducted to evaluate the precision and performance of these models in predicting geoid heights.The comparison of free-air gravity anomalies, geoid undulation values, and residual geoid heights among the different GGMs reveals different levels of accuracy and performance.XGM2019e emerges as the most accurate model in capturing both free-air gravity anomalies and geoid undulation values, exhibiting the lowest RMSE values and relatively small maximum differences compared to actual measurements.Tongji-Grace02k also demonstrates promising results, with low mean values and RMSE values, despite having slightly lower mean residuals.However, Tongji-GMMG2021S exhibits the highest RMSE and standard deviation values, indicating larger deviations from actual measurements.Overall, this comparison highlights the importance of selecting the appropriate GGM for specific applications and emphasizes the ongoing requirement to refine and enhance the GGMs to improve their accuracy and precision of prediction.The study concludes that the geoid model derived from XGM2019e offers superior accuracy, with standard deviations of residuals at 0.306 meters, a geoid undulation RMSE of ±0.033 meters, and the minimum free-air gravity anomalies at 10.291 mGal, compared to another model.would like to express our gratitude to the DSMM for supplying BM/SBM values, and airborne and terrestrial gravity data for the Peninsular Malaysia.Furthermore, we would like to extend our appreciation to ICGEM, NASA, and DTU Space National Space Institute for their generous provision of the respective global geopotential model, SRTM, and DTU21 ocean-wide gravity field datasets.

Figure 2 .
Figure 2. Distribution of (a) terrestrial gravity points, (b) airborne gravity points, (c) marine gravity points, and (d) Global DEM from the SRTM model.

Figure 5 .
Figure 5. Development of gravimetric geoid model by using KTH method.