Hong Kong Height Values Transformation

The 2D Method - 7-parameters trnasformation

GIS projects in Hong Kong often require the use of GNSS to collect your location and report the results in the Hong Kong 1980 Grid (EPSG: 2326). If Python is your tool, PyProj is the simplest way to transform your location between latitude, longitude to HK1980 Grid Northing and Easting. It has been well documented and is very easy to use.

The method used to transform coordinates is known here as the 7-parameters transformation, or the Helmert transformation. It rotates, scales and translates one set of coordinates to a different projection.

from pyproj import Transformer
from pyproj import CRS

epsg4326 = CRS.from_epsg(4326) # WGS1984
epsg2326 = CRS.from_epsg(2326) # Hong Kong 1980 Grid

x = 114.215286
y = 22.285663
z = 10

transformer = Transformer.from_crs(epsg4326, epsg2326)
# deliberate swapping of y and x
nx,ny,nz=transformer.transform(y,x,z) 

print("EPSG2326", nx,ny,nz)

However, this method does not transform the height values. See this document for more details Lands Department Hong Kong - Geodetic Datum Transformation.

Using the Web-based Transformation Tool provided by the Lands Department of Hong Kong, we could test out some transformation results.

Scratching the Surface ~wink~

The ellipsoid is a representation of the earth's surface, or generalized height model of the earth as a whole. Ellipsoidal Height is the height above this ellipsoid and the height value measured from your GNSS.

The geoid is a better representation of local surface height values or a region, as it is more detailed than the ellipsoid. The Orthometric Height is the distance between the geoid height and the measured point. The orthometric height that we adapt in Hong Kong is the height above Hong Kong Principal Datum (HKPD).

There is plenty of well written articles on these confusing surfaces and heights. Here is one, if you would like to know more Geoid vs. Ellipsoid: What’s the Difference?.

Ellipsoidal height to HKPD height

The Lands Department does provide a list of control points height values. The source can be obtained here. These control points are used to calculate the height model of Hong Kong, helping us match the height value of ellipsoidal height to height above HKPD. Don't ask me why they published a PDF instead of a machine readable format, but the good news is there are only a few dozen points.

Below is a snippet of the CSV that I've collected (separated by space), through copy and pasting the PDF file above. Refer to my repository for the full file.

control_pts_heights.csv
141 22 14 59.09628 114 10 28.81741 11.921 14.460
142 22 15 21.71561 114 12 44.85585 152.023 154.378
143 22 11 56.64407 114 13 13.50645 159.095 161.318
146 22 20 44.86275 114 12 54.42899 268.501 271.065
...

The calculations steps are simple:

1. Read the CSV and find the height difference for each control point
2. Find the 4 closest control point stations and mean the height difference
3. Add the mean height difference to the input height value

And here is the code to run the calculations:

# read file using pandas without header and convert it to dataframes (control_pts_heights.csv)
names=["stn", "lat_d", "lat_m", "lat_s", "lon_d", "lon_m", "lon_s", "ell_h", "pd_h"]
csv = pd.read_csv('control_pts_heights.csv', header=None, sep=" ").values
df = pd.DataFrame(csv, columns=names)

# calculate each station's height difference (hd) between the two heights
df['lat'] = df.apply(dms2lat, axis=1)
df['lon'] = df.apply(dms2lon, axis=1)
gdf = gpd.GeoDataFrame(df, geometry=gpd.points_from_xy(df["lon"], df["lat"]))
gdf["hd"] =  gdf["pd_h"] - gdf["ell_h"]

# define lat lon height
lat = 22.285663
lon = 114.215286
height = 10.0

# finding the 4 closest stations and add the mean height difference the input height value
stations_kd_map = point_to_closest_stations_id((lon, lat))
stations_keys = [k for k,v in stations_kd_map]
stns = gdf.query('stn in @stations_keys')

n_height = height + stns["hd"].mean()
print(n_height) # Northing and Easting comes from the step before (7 parameters)

2-Step Approach - Combining the Two

The full transformation procedure requires two steps:

1. 7-parameters transformation of latitude and longitude to Northing and Easting
2. Ellipsoidal height to HKPD height

Both steps are required to obtain an accurate 3D position transformation. You may reverse the steps to transform from Hong Kong 1980 Grid to WGS84.

Lastly...

Since not everyone has a professional software available (the local government only provided a guide for Leica LGO), I hope this article could help another fellow mapper and reduce your development time.

Please refer to my repository to run the transformation along with the prepared station CSV.

Further Reading

The details of all the theories and detailed guide using Leica LGO Version 3.0 to process this height projection change can be found in this link: Suggested Procedure of Using Leica LGO Software for Transforming GPS Surveyed Height (Ellipsoidal Height) to Local Height (Height above HKPD) within a Small Local Area.