The worldwide physical height datum project

The definition of a common global vertical coordinate system is nowadays one of the key points in Geodesy. With the advent of GNSS, a coherent global height has been made available to users. The ellipsoidal height can be obtained with respect to a given geocentric ellipsoid in a fast and precise way using GNSS techniques. On the other hand, the traditional orthometric height is not coherent at global scale. Spirit levelling allows the estimation of height increments so that orthometric heights of surveyed points can be obtained starting from a benchmark of known orthometric heights. As it is well known, this vertical coordinate refers to the geoid, which is assumed to be coincident to the mean sea level. By means of a tide gauge, the mean sea level is estimated and thus a point of known orthometric height is defined. This assumption, which was acceptable in the past, became obsolete given the level of precision which is now required. Based on the altimetry observation, one can precisely quantify the existing discrepancy between geoid and mean sea level that can amount to 1 ÷ 2 m at global scale. Therefore, different tide gauges provide biased estimates of the geoid, given the discrepancy between this equipotential surface and the mean sea level. Also, in the last years, another vertical coordinate was used, the normal height that was introduced in the context of the Molodensky theory. In this paper, a review of the existing different height systems is given and the relationships among them are revised. Furthermore, an approach for unifying normal height referring to different tide gauges is presented and applied to the Italian test case. Finally, a method for defining a physical height system that is globally coherent is discussed in the context of the definition of the International Height Reference System/Frame, a project supported by the Global Geodetic Observing System of the International Association of Geodesy (IAG). This project was established in 2015 during the XXVI IAG General Assembly in Prague as described in IAG Resolution no. 1 that was presented and adopted there.


Introduction
The definition and the realization of a physical vertical datum and its connection to the mean sea level (MSL) have been central problems in Geodesy. The discussion on these relevant topics traces back to the seventies (see, e.g. Balazs 1973 ;Sturges 1974;Lelgemann 1977), with further developments in the 80s (see, e.g. Colombo 1980;Rummel and Teunissen 1988). In more recent years, several works have analysed the physical height datum unification problem. The methods proposed for datum unification exploit the availability of precise geoid models and GNSS observations at levelling benchmarks (see, e.g. Ardalan et al, 2002;Amos and Featherstone 2009) or are discussed in the theoretical frame of the Geodetic Boundary Value Problem (see, e.g. Xu 1992;Nahavandchi and Sjoberg 1998;Sansò and Venuti 2002;Sacerdote and Sansò 2004;Sjoberg 2011). The availability of satellite only Global Geopotential Models (GGM) has then opened new perspectives on the height datum unification problem as described, among the others, in Rummel (2012), Gatti et al. (2013) and Amjadiparvar et al. (2013) where the contribution of the GOCE mission (Drinkwater et al. 2003) is analysed.
This paper deals with the definition of a global physical height datum. In Sect. 2, the different definitions of physical 1 3 and geometrical heights are reviewed and discussed. In Sect. 3, the method for solving the physical height unification problem described in Barzaghi et al. (2015) is revised. In Sect. 4, a possible approach for establishing a global physical height system is illustrated. This is done in the framework of the activities of Focus Area "Unified Height System" of the Global Geodetic Observing System (GGOS) of the International Association of Geodesy (IAG) and of the International Height Reference System/Frame (IHRS/IHRF) project (Ihde et al. 2017;Sànchez 2013;Sànchez and Sideris 2017). Conclusions are then given in Sect. 5.

The different definitions of height
Basically, we can define two kinds of heights, i.e. the physical heights defined in terms of the Earth gravity potential W(P) and the geometric heights from GNSS that are only weakly dependent on the gravity field (via the satellite orbits). As said, the physical height definition is strictly connected to the gravity potential of the Earth W(P) . In a given point, W(P) is the sum of the gravitational potential V(P) due to the Earth masses and of the centrifugal potential CP(P) due to the Earth rotation. Thus, we can write where (Q) is the Earth mass density in the Earth volume , is the angular velocity of the Earth and x P , y P are the cartesian coordinates of P in a geocentric reference frame. The definition of the physical heights is based on that of geopotential number C P where W 0 = W(P 0 ) , with P 0 on the geoid, and W P = W(P), with P on the Earth surface. Since, along the plumb line, it holds that with g = |∇W| , we can also write where the integral is taken along the plumb line from P 0 to P. (1) By dividing C P for a conventional value 0 = |∇U| (e.g. with Finally, the normal height can be defined. For that, we consider the normal potential U and the modulus of the normal gravity = |∇U| . If we assume that on the ellipsoid U 0 = W 0 and in a point Q along the plumb line U Q = W P , we have where the integral is along the ellipsoidal normal. The normal height H * is then defined as with Thus, we can say that the normal height of P , H * (P) , is the height of Q above the ellipsoid. The definitions given above show that the different height coordinates are strictly related in the sense that they all depend on C P , i.e. on the gravity potential. Among them, the most frequently used are H and H * . Particularly, H * are used in Europe for the definition of the European Vertical Reference System (EVRS) (Ihde and Augath 2001).
H dyn , H and H * can be obtained by spirit levelling increments properly corrected for gravity. As it holds that One can get for H dyn the so-called dynamic correction ΔH dyn : with ΔL AB = ∫ B A l measured levelling increment and ΔH dyn dynamic correction.
We can also see that, by (12), geopotential number differences are obtained from levelling increments as Similarly, orthometric correction ΔH and normal corrections ΔH * can be obtained. The orthometric correction is derived starting from the following identity: where 0 is a given constant value of the normal gravity. Hence, we have, taking (14) into account where ΔL AB is, as in (12), the measured levelling increments and ΔH ort is the orthometric correction. In the same way, we can get normal height increments from ΔL AB as where ΔH * is the normal correction.
Furthermore, it has to be mentioned that, in recent years, due to the availability of atomic clocks, C P (and thus heights) can be determined via the general relativity formula (Shen et al. 2011) (12) In (17), f P is the oscillator frequency at P having potential W P while f 0 is its frequency in P 0 with potential W 0 ; c is the speed of light. Thus, by highly precise observations of the frequency shift, C P can be obtained.
Besides the physical heights discussed above, geometric heights are nowadays frequently used. As stated at the beginning of this section, geometric heights are obtained from GNSS observations. The geometric height of a point P is the ellipsoidal height h(P) . It refers to the ellipsoid and is defined as the distance from the point P on the surface and the point P 0 onto the ellipsoid along the perpendicular to the ellipsoid in P . As already pointed out, although it is claimed that is a purely geometric quantity, a weak connection to gravity is present. This is through the satellite orbits that enter in the estimation of the point coordinates obtained via GNSS techniques. The geometric heights h(P) are then linked to the physical heights H(P) and H * (P) via the very well known formulas (Heiskanen and Moritz, 1967;Sansò et al., 2019) where N(P) is the geoid undulation and (P) is height anomaly.

Physical height system unification: estimating the biases
The traditional definition of a national height system is based on spirit levelling starting from a tide gauge. The assumption is that the tide gauge benchmark refers to the geoid. As it is well known, this is not a correct statement. The tide gauge benchmark has a known height with respect to the mean sea level, which is close to the geoid but not coincident with it. Since, as said, the discrepancy between geoid and mean sea level can reach 1 ÷ 2 m, inconsistencies among physical height systems referred to different tide gauges can have the same order of magnitude. This has been accepted in the past since the global unification of the physical height systems faced some practical difficulties. Merging different national height systems requires a common adjustment of national levelling networks which are to be referred to a unique tide gauge. Although in recent years this has been done at continental level (e.g. in Europe with the EVRS, Ihde and Augath 2001), this cannot be managed at global level. Due to the availability of GNSS techniques, the unification problem can be tackled in a different way (Gatti et al. 2013;Barzaghi et al. 2015). By Eqs. (8) and (9), in a given tide gauge point P j 0 to which heights in the j th area refer, we can write that where W j is the difference between the gravity potential in P j 0 and its value W 0 on the geoid. If we now consider the ellipsoidal heights from GNSS, we can write in any given point P of the j th area the following equation: as, by definition, h(P) − H * (P) = (P) and by Bruns' formula (P) = T(P) (P) (Heiskanen and Moritz 1967). Then by introducing the biased height anomaly However, an important remark on the T(P) values has to be stressed. T(P) cannot be estimated using local gravity data because they are biased too since height information are used in observing and reducing them. For the same reason, high d/o GGMs are not suitable. As a matter of fact, ground gravity data are used in estimating these models and thus the problem detailed above holds also for them. One possible solution of this problem is to use a combination of a satellite only GGM and of a high d/o GGM. In this way, T(P) can be represented as In (24), T L (P) could be given as one of the GOCE models to d/o 200 and T H (P) as a high-degree model, e.g. EGM2008 from d/o 201 to full resolution (Pavlis et al. 2012). It has been proved that in this way the impact of the bias effect in T(P) is minimized (Gatti et al. 2013). This can be heuristically explained in the following way. The T L (P) is bias independent since no ground data have been used for estimating it. On the other hand, the high-frequency component T H (P) is biased but has a limited impact on T(P) since this high-frequency component is less than 10% of the whole T(P) signal. Based on the assumption in (24), we can then write the observation, Eq. (23), as where ∼ j (P) is an observed value coming from combined GNSS and spirit levelling measurements. The covariance structure of the observed value in the left hand of (25) must be also carefully detailed. By assuming independence of the three terms, we have The covariance of ∼ j (P) can be set equal to with 2 ∼ obtained by error propagation from GNSS and levelling observations. C T L can be computed as in, e.g. Gerlach and Fecher (2012) propagating the block diagonal structure of the GOCO model. Finally, C T H can be estimated from, e.g. EGM2008 global model error geographically rescaled (Gatti et al. 2013).
A numerical test of this procedure has been performed in Italy. In fact, the levelling lines of the peninsular part of Italy refer to the main tide gauge in Genova while those of Sicily and Sardinia refer to tide gauges in Catania and Cagliari, respectively. Thus, the Italian levelling network does not refer to a common reference point and the proposed unification method based on Eq. (25) can be applied. The results of this test are described in Barzaghi et al. (2015). The estimated biases, of the order of 10 ÷ 20 cm, are coherent with the magnitude of the Mean Dynamic Ocean Topography of the Mediterranean Sea (Barzaghi et al. 2009). Furthermore, the bias between Italy mainland and Sicily is consistent with a direct estimate based on spirit and trigonometric levelling obtained by the Istituto Geografico Militare. Some other approaches similar to the one described in this section are currently used. They are based on the same hypotheses assumed here and are documented in, e.g. Fotopoulos (2013).

The global physical height system definition according to the IHRS/IHRF approach
This approach for unifying the physical heights has been proposed during the IAG/IUGG assembly held in 2015 in Prague, Czech Republic. It has been adopted by the International Association of Geodesy (IAG) as Resolution 1: Definition and Realization of an International Height Reference System (IHRS).
In that resolution, it was stated to establish a physical height system following these conventions: (1) The vertical reference level is an equipotential surface of the Earth gravity field with potential value W 0 (at the geoid).
(2) Parameters, observations and data shall be related to the mean tidal system/mean crust. (3) The unit of height is the metre and the unit of time is the second (IS). (4) The vertical coordinates are the differences −ΔW P between the potential W P of the Earth gravity field at the considered point P and the geoidal geopotential value W 0 ; the potential difference −ΔW P is also designated as C P : −ΔW P = C P = W 0 − W P . (5) The spatial reference of the position P for the potential It is then further stated that the conventional value of W 0 is More details on IHRS can be found in Ihde et al. (2017). The realization of the IHRS is now in progress. In the years, the International Height Reference Frame (IHRF), i.e. the realization of IHRS, will be established as a global network of points with known C P values. The aim is to define C P at 1 cm level in terms of geoid undulation/height anomaly. This process is nowadays carried out by the Global Geodetic (28) W 0 = 62636853.4m 2 s 2 .
Observing System, the International Gravity Field Service of IAG and IAG Commission 2. The first design of the IHRF network has been proposed at the IAG/IUGG General Assembly in Montreal, Canada, in 2019 and is described in Sànchez et al. (2020).
With this approach, it is clearly recognized that any physical height system is basically derived as a function of W(P), as is stated in Sect. 2. Thus, to establish the IHRS, one has to give a proper estimate of W(P) at any given point P . Usually, the gravity potential in P is defined as where U(P) is the normal potential of the reference ellipsoid GRS80 and T(P) is the anomalous potential (Heiskanen and Moritz 1967). Since in IAG Resolution 1 the W 0 reference value is different from U onto the GRS80 ellipsoid, U GRS80 o , one has to consider the modified equation: where the ΔW 0 term accounts for this discrepancy. Given that U(P) can be analytically computed (Moritz 1980) and that U 0 and W 0 are known, the estimation of W(P) is equivalent to the estimation of T(P).
In Geodesy, there are very well known methods for estimating T(P) (Sansò and Sideris 2013;Sjöberg 2003;Vanicek et al. 2013). A straightforward way of getting T(P) is to use a high d/o GGM. However, in view of what we discussed in Sect. 3, this value is biased since high d/o GGMs contain ground (biased) gravity data. As a matter of fact, some tests performed in Sànchez et al. (2020) proved that biases at the level of decimeters in terms of geoid undulation exist among T(P) estimates from different GGMs in a given point.
So, to properly estimate T(P) , we should rely on satellite only GGMs and newly collected local gravity data having GNSS surveyed coordinates. Following this approach, one possible scheme for estimating T(P) at a point P on the Earth surface can be designed as This is the standard remove-compute-restore adopted in the Molodensky context and based on collocation for the estimation of T(P) . In this scheme, no biased height information enters. The g(P) values are not biased since GNSS heights h(P) are not affected by the tide gauge problem described in Sect. 3. The same holds for the satellite only GGM that does not contain ground gravity information. Similarly, the SRTM derived DTM that is frequently used in the RTE computation can be given in terms of ellipsoidal height. Thus, the g RTE (P) is bias free too.
Finally, T (P) is estimated by applying Bruns' equation, which depends on ̂ (P) and the known (P) value computed at the ellipsoidal height h(P) . Hence, by a proper combination of gravity data, satellite data and GNSS positioning, we can obtain the unbiased estimation of T(P) and, by (30), of W(P).

Conclusions
The physical height datum problem has been thoroughly studied from the theoretical point of view and different approaches have been devised to solve it. However, in recent years, GNSS techniques allowed the estimate of the ellipsoidal heights in a fast and reliable way and the problem of getting physical heights on a global scale was considered less important in the geodetic community.
Nevertheless, physical heights are still important since they are strictly related to the gravity potential and are suitable when dealing with problems like water flows. As an example, the problem of estimating the mean dynamic ocean topography (and its changes in time) at global scale can be mentioned.
In the last decades, mainly due to availability of a new generation of satellite only GGMs, the problem of estimating a global physical height system was revisited and new efforts were performed to practically implement it. In this paper, two possible methods for physical height system unification are presented. The first method has been widely applied and is based on the estimation of the bias between the gravity potential W 0 at the geoid and the gravity potential in a point P on the mean sea level as measured at a given tide gauge. This bias is different from one tide gauge to another so that heights referred to different tide gauges are not coherent. The solution proposed in the paper allows the estimate of these biases based on the availability of spirit levelling, GNSS data and a proper GGM. An application of this method to the Italian case proved its feasibility with an estimate of the bias between the physical heights in the peninsular part of Italy and those in Sicily close to the direct estimate obtained by IGM. Even though this method is suitable for the unification of different physical height systems at global scale, a more general solution has been recently proposed in the frame of the IHRS/IHRF project of IAG. This approach is based on the estimate of the gravity potential of the Earth W(P) in any given point P and the assumption of a conventional value of W 0 , the gravity potential on the geoid. In this way, the geopotential numbers C P can be computed and orthometric and normal heights can be derived. The scheme proposed in this paper, based on a satellite only GGM, local gravity data and GNSS observations, can lead to the unbiased estimate of the Earth gravity potential W(P) with a precision, in terms of geoid/height anomaly, of the order of 2 ÷ 3 cm at global scale. Thus, in this way, a global physical height system can be obtained that improves the present day situation by two orders of magnitude.
Funding This research received no external funding. Open access funding provided by Politecnico di Milano within the CRUI-CARE Agreement.

Compliance with ethical standards
Conflict of interest The authors declare no conflict of interest.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creat iveco mmons .org/licen ses/by/4.0/.