作者:郜科科

两个坐标系统的参考椭球不同,实地一个点的不同坐标系的值是不同的,不同的部门采用的坐标系统经常是不一致,所以要转换后才能相互利用。例如目前使用的北京市观测站点位置根据GPS的定位而来,GPS使用的地理坐标系为GCS_WGS_1984,所以其坐标的地理坐标系也为GCS_WGS_1984,而假如需要将这些点显示在Web端的地图上,Web端的投影坐标系WGS_1984_Web_Mercator_Auxiliary_Sphere,就需要进行地理坐标系转换为投影坐标系的操作。

关于地理坐标系投影坐标系的关系这里不再赘述,这里介绍一下WKID。WKID的英文全称是Well Known ID,即众所周知的编号。这个编号是大家坐下来一起讨论、约定和认同的,具体有唯一性。众多的坐标系统有了自己的WKID,就像每个人都有自己的身份证号一样,从出生就定了,即使是名字改了,还是可 能通过身份证号确定,这为空间数据的使用、转换、共享等起到关键作用。

例如下表为GCS_WGS_1984的WKID格式及某点的投影文件的具体实例:

WKID
4326
名称
GCS_WGS_1984
参数
GEOGCS["GCS_WGS_1984",DATUM
["D_WGS_1984",SPHEROID["WGS_1984",6378137.0,298.257223563]],
PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]]
需要转换为的投影坐标同样有其WKID,例如下边为WGS_1984_Web_Mercator_Auxiliary_Sphere的具体实例:
WKID
3857
名称
WGS_1984_Web_Mercator_Auxiliary_Sphere
参数
PROJCS["WGS_1984_Web_Mercator_Auxiliary_Sphere",GEOGCS
["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137.0,298.257223563]],
PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],PROJECTION["Mercator_Auxiliary_Sphere"],
PARAMETER["False_Easting",0.0],PARAMETER["False_Northing",0.0],PARAMETER["Central_Meridian",0.0],
PARAMETER["Standard_Parallel_1",0.0],PARAMETER["Auxiliary_Sphere_Type",0.0],UNIT["Meter",1.0]]
如果对自己需要进行转换的坐标系的WKID不了解,可以从以下两个网站进行查询:
下面列举一些我国常用的地理坐标系的WKID:
坐标系
WKID
备注
地理坐标
4214
GCS_Beijing_1954
地理坐标
4326
GCS_WGS_1984
地理坐标
4490
GCS_China_Geodetic_Coordinate_System_2000
地理坐标
4610
GCS_Xian_1980

进行坐标系的转换有很多工具,其中比较常用的又Proj.4相关库,如果我们使用Python进坐标转换的话,有高级的Pyproj第三方库可以使用。其文档地址如下:

这个库非常简单,我们只需要掌握其中的一个主要函数就可以了:

transform(p1, p2, x, y, z=None, radians=False)
示例:x2, y2, z2 =transform(p1, p2, x1, y1, z1, radians=False)

这个函数表示在p1坐标系和p2坐标系之间进行坐标转换,x1,y1,z1是由p1坐标系定义的坐标,z为高度单位是米。X2,y2,z3是由p2坐标系定义的坐标,它是经过转换过后返回的,默认z1=none。Radians参数表示是否用弧度返回值。

下面我们进行一下北京市观测站点的坐标转换,如下所示为转换的代码:

# 投影变换
def proj_trans():
# 读取经纬度
data = pd.read_excel(u"D:/Visualization/python/file/location.xlsx")
lon = data.lon.values
lat = data.lat.values
print lon, lat
p1 = pyproj.Proj(init="epsg:4326")  # 定义数据地理坐标系
p2 = pyproj.Proj(init="epsg:3857")  # 定义转换投影坐标系
x1, y1 = p1(lon, lat)
x2, y2 = pyproj.transform(p1, p2, x1, y1, radians=True)
print x2, y2

在上述代码中需要注意需要在转换前首先定义数据的地理坐标系和转换后的投影坐标系,这样才能进行有目的性的转换。

转换前后的结果如下所示:

name
lon
lat
x
y
万寿西宫
116.366
39.8673
12953803.87
4846677.374
定陵
116.17
40.2865
12931985.25
4907663.441
东四
116.434
39.9522
12961373.59
4858998.543
天坛
116.434
39.8745
12961373.59
4847721.686
农展馆
116.473
39.9716
12965715.05
4861816.127
官园
116.361
39.9425
12953247.27
4857590.051
万柳
116.315
39.9934
12948126.57
4864983.232
顺义
116.72
40.1438
12993210.97
4886860.96
怀柔
116.644
40.3937
12984750.68
4923319.714
昌平
116.23
40.1952
12938664.41
4894348.889
奥体中心
116.407
40.0031
12958367.96
4866392.773
古城
116.225
39.9279
12938107.82
4855470.429

下面我们使用ArcMap中的Project工具进行实验,并对照一下实验结果,如下所示为在ArcMap中进行投影转换后的结果,与上述结果基本相同。


标签:1984,转换,python,WKID,Pyproj,GCS,坐标,WGS,坐标系