GeoPandas to_crs()转换坐标时Z值丢失如何保留POINT Z类型

基于存放在metadata_df变量中的影像元数据CSV，提取经纬度、高程字段构建Point几何生成GeoDataFrame的初始代码可正常运行：

df = metadata_df
epsg_code='4979'

geometry = [Point(xyz) for xyz in zip(df['lon'], df['lat'], df['alt'])]     
gdf1 = gpd.GeoDataFrame(df, geometry=geometry, crs='EPSG:'+epsg_code)

生成的gdf1坐标系为EPSG:4979，几何为带Z值的POINT Z类型，前5条几何打印结果正常：

print(gdf1.geometry.head())
0    POINT Z (-110.72734 43.71498 4629.18463)
1    POINT Z (-110.72974 43.74431 4625.70121)
2    POINT Z (-110.73011 43.77137 4614.43130)
3    POINT Z (-110.73102 43.80133 4612.92112)
4    POINT Z (-110.73159 43.83043 4598.42577)
Name: geometry, dtype: geometry

调用to_crs()方法将其转换为EPSG:4978坐标系时Z值被完全剥离，几何退化为二维POINT类型：

gdf2 = gdf1.to_crs("EPSG:4978")
print(gdf2.geometry.head())
0    POINT (-1634202.022 -4318557.192)
1    POINT (-1633585.179 -4316381.716)
2    POINT (-1632877.279 -4314426.871)
3    POINT (-1632130.428 -4312246.197)
4    POINT (-1631381.075 -4310136.444)
Name: geometry, dtype: geometry

当前运行环境依赖版本为：

GeoPandas 0.9.0
PyProj 3.3.1
PyGeos 0.12

需要实现EPSG:4978坐标转换过程中保留Z值，让to_crs()输出POINT Z类型而非二维POINT。

GeoPandas 0.9.0版本对三维几何的重投影支持存在缺陷，to_crs()方法默认仅处理X、Y两个维度的坐标转换，会直接丢弃Z维度值，哪怕源坐标系和目标坐标系均为带高程定义的三维坐标系（EPSG:4979为WGS84三维大地坐标系，EPSG:4978为WGS84地心地固三维坐标系，原生支持三维坐标定义）。

可根据环境限制二选一：

• 方案一（推荐）：升级GeoPandas版本至0.10.0及以上。高版本GeoPandas已原生支持三维几何的全维度坐标重投影，无需修改原有业务代码，直接调用to_crs("EPSG:4978")即可自动保留Z值，输出正确的POINT Z类型结果。
• 方案二（无法升级GeoPandas时使用）：绕过GeoPandas内置的to_crs()方法，直接调用PyProj的三维坐标转换接口手动处理几何，实现代码如下：

from pyproj import Transformer
from shapely.geometry import Point

# 初始化三维坐标转换器，always_xy=True指定输入坐标顺序为经度、纬度、高程，避免轴顺序错误
transformer = Transformer.from_crs("EPSG:4979", "EPSG:4978", always_xy=True)

def convert_3d_point(point):
    x, y, z = transformer.transform(point.x, point.y, point.z)
    return Point(x, y, z)

# 批量转换几何并赋值坐标系
gdf2 = gdf1.copy()
gdf2["geometry"] = gdf1.geometry.apply(convert_3d_point)
gdf2.crs = "EPSG:4978"

转换完成后打印gdf2.geometry.head()即可验证，输出结果为带正确Z值的POINT Z类型，无维度丢失问题。

内容的提问来源于stack exchange，提问作者E. Case