在GIS基础算法库中,有很多工具都支持投影变化。正是使用投影变化算法,实现不同坐标系之间数据的转换。
确认坐标系
在执行投影变化前,首先需要确认你的投影坐标使用的坐标系。
如果你对投影坐标数据使用的坐标系不太清楚的话,请参考这篇文章,尝试猜测数据坐标系。
pyproj单点转换
本文选择pyproj库提供的Transformer类。
例如,投影坐标为x=-11705384.6174 y=4826113.6622。投影坐标系为EPSG:3857。经纬度坐标使用EPSG:4326坐标系。示例代码,如下所示:
from pyproj import Transformer
transformer = Transformer.from_crs("EPSG:3857", "EPSG:4326")
latitude, longitude = transformer.transform(-11705384.6174, 4826113.6622)
注意:以上代码得到的经纬度顺序与投影坐标XY不一致,进行了交换。如果不希望交互经纬度的前后顺序,请使用参数always_xy=True。代码如下所示:
transformer = Transformer.from_crs("EPSG:3857", "EPSG:4326", always_xy=True)
longitude, latitude = transformer.transform(-11705384.6174, 4826113.6622)
pyproj多点批量转换
如果需要批量转换,参考如下示例:
import numpy as np
# 使用numpy数组进行高效批量转换
x_array = np.array([-11705384.6174, -11705385.0, -11705386.0])
y_array = np.array([4826113.6622, 4826114.0, 4826115.0])
latitudes, longitudes = transformer.transform(x_array, y_array)
pyproj自定义坐标系
在定义坐标系时,pyproj除了支持EPSG方式,还支持proj字符串和WKT格式,方便构建自定义投影坐标系。
proj字符串定义示例:
proj_wgs84 = "+proj=longlat +datum=WGS84 +no_defs"
proj_utm50 = "+proj=utm +zone=50 +datum=WGS84 +units=m +no_defs"
transformer = Transformer.from_crs(proj_utm50, proj_wgs84)
WKT定义示例:
wkt_3857 = 'PROJCS["WGS 84 / Pseudo-Mercator", GEOGCS["WGS 84", DATUM["WGS_1984", SPHEROID["WGS 84", 6378137, 298.257223563]], PRIMEM["Greenwich", 0], UNIT["degree", 0.0174532925199433]], PROJECTION["Mercator_1SP"], PARAMETER["central_meridian", 0], PARAMETER["scale_factor", 1], PARAMETER["false_easting", 0], PARAMETER["false_northing", 0], UNIT["metre", 1]]'
transformer = Transformer.from_crs(wkt_3857, "EPSG:4326")
更多pyproj使用问题,欢迎留言或联系我们。转载须注明出处。