13

如何使用Python读取GeoJSON和Shapefile文件，并转换坐标系

 3 years ago
source link: https://zhuanlan.zhihu.com/p/356822042
Go to the source link to view the article. You can view the picture content, updated content and better typesetting reading experience. If the link is broken, please click the button below to view the snapshot at that time.

如何使用Python读取GeoJSON和Shapefile文件，并转换坐标系

上海市新能源汽车公共数据采集与监测研究中心 数据分析师

下面这个网页可以不写代码进行Shapefile和GeoJSON的坐标系转换

Shapefile是美国环境系统研究所公司（ESRI）开发的一种空间数据开放格式，而GeoJSON是一种基于JSON对地理数据结构进行编码的格式。Shapefile和GeoJSON是目前使用最频繁的两种地理数据格式。

由于国内存在多种不同的地理坐标系（如WGS-84、GCJ-02、BAIDU-09等），因此在使用这两种格式的地理数据时经常需要进行坐标系转换的操作。

本文讲述如何使用Python读取Shapefile和GeoJSON文件，并转换地理坐标系。

以下图所示的的上海建筑边界数据为例，格式为Shapefile，坐标系为WGS-84，底图为高德地图，坐标系是GCJ-02，可以看到Shapefile数据与地图是不匹配的。

v2-eb68c085371f926df03ab207ab64d651_720w.jpg

第一步：使用Python中的Geopandas插件读取Shapefile文件

# 导入Geopandas
import geopandas as gp

# 读取Shapefile文件
loadData = gp.read_file("E:\\地图数据归档\\AOI数据\AOI_GCJ02\\aoi_gcj.shp")
print("读取的数据：")
print(loadData)

# 判断Shapefile文件的格式，有点、线、面三种
loadType = loadData['geometry'][0].geom_type
print("空间数据类型：")
print(loadType)

第二步：读取空间数据中的每一个点，并进行地理坐标系转换

不同类型的空间数据有不同的存储结果，因此需要用不同的代码进行坐标系转换

# 引入地理坐标系转换功能包
from gc_method import CoorChange as cc
# 引入Point类型
from shapely.geometry import Point

originCoordSystem = "GCJ-02"  # 转换前地理坐标系
afterCoordSystem = "WGS-84"  # 转换后地理坐标系

# 循环读取每个Geometry
for i in range(0, len(loadData)):
    # 读取当前点数据
    currentGeometry = loadData.loc[i, 'geometry']
    # 读取当前点的经度和纬度
    currentLon = currentGeometry.x
    currentLat = currentGeometry.y
    # 进行地理坐标系转换
    transformedResult = cc.main(currentLon, currentLat, originCoordSystem, afterCoordSystem)
    # 生成转换之后的点数据
    transformedPoint = Point(transformedResult)
    # 用转换后的坐标替换原来的坐标
    loadData.loc[i, 'geometry'] = transformedPoint
# 引入地理坐标系转换功能包
from gc_method import CoorChange as cc
# 引入Point类型
from shapely.geometry import Point,LineString

originCoordSystem = "GCJ-02"  # 转换前地理坐标系
afterCoordSystem = "WGS-84"  # 转换后地理坐标系

# 循环读取每个Geometry
for i in range(0, len(loadData)):
    # 读取当前线数据
    currentGeometry = loadData.loc[i, 'geometry']
    # 读取线数据中的每个坐标点
    currentLine = currentGeometry.coords
    transformedList = []
    for j in range(0, len(currentLine)):
        # 进行坐标系转换
        transformedResult = cc.main(currentLine[j][0], currentLine[j][1], originCoordSystem, afterCoordSystem)
        transformedList.append(Point(transformedResult[0], transformedResult[1]))
    # 生成转换之后的线数据
    transformedLine = LineString(transformedList)
    # 用转换后的线数据替换原来的线数据
    loadData.loc[i, 'geometry'] = transformedLine
# 引入地理坐标系转换功能包
from gc_method import CoorChange as cc
# 引入Point类型
from shapely.geometry import Point,Polygon

originCoordSystem = "GCJ-02"  # 转换前地理坐标系
afterCoordSystem = "WGS-84"  # 转换后地理坐标系

# 循环读取每个Geometry
for i in range(0, len(loadData)):
    # 读取当前面数据
    currentGeometry = loadData.loc[i, 'geometry']
    # 读取当前面的每个坐标点
    currentPolygon = currentGeometry.exterior.coords
    transformedList = []
    for j in range(0, len(currentPolygon)):
        # 进行坐标系转换
        transformedResult = cc.main(currentPolygon[j][0], currentPolygon[j][1], originCoordSystem,
                              afterCoordSystem)
        transformedList.append(Point(transformedResult[0], transformedResult[1]))
    # 生成转换之后的面数据
    transformedPolygon = Polygon(transformedList)
    # 用转换后的面数据替换原来的面数据
    loadData.loc[i, 'geometry'] = transformedPolygon

地理坐标系转换的具体代码可参考下面这篇文章：

第三步：将坐标系转换后的数据以Shapefile或GeoJSON格式保存

# 输出的数据格式
fileType = 'shapefile'

# 进行保存
if fileType == 'geojson':
    loadData.to_file("E:\\地图数据归档\\AOI数据\AOI_GCJ02\\aoi_gcj_after.geojson", driver='GeoJSON', encoding="GB2312")
elif fileType == 'shapefile':
    loadData.to_file("E:\\地图数据归档\\AOI数据\AOI_GCJ02\\aoi_gcj_after.shp", driver='ESRI Shapefile',
                      encoding='GB2312')

下图为坐标系转换后的数据

下面这个网页可以不写代码完成GeoJSON和Shapefile的地理坐标系转换

v2-2103250f5c47784db85fc2f95193c451_720w.jpg

如果觉得文本对你有帮助，请点赞、收藏并关注，谢谢！


report

Recommend

About Joyk


Aggregate valuable and interesting links.
Joyk means Joy of geeK