跳到主要内容

投影

SpatialReference​

操作投影​

创建投影​
# 创建投影
prosrs: osr.SpatialReference = osr.SpatialReference()
prosrs.ImportFromWkt(prj)
prosrs.ImportFromEPSG(4326)
# 获取投影对应地理坐标系
geosrs = prosrs.CloneGeogCS()
设置投影策略​
# 设置坐标策略
# osr.OAMS_TRADITIONAL_GIS_ORDER, 即 X 对应 lng, Y 对应 lat
prosrs.SetAxisMappingStrategy(osr.OAMS_TRADITIONAL_GIS_ORDER)
geosrs.SetAxisMappingStrategy(osr.OAMS_TRADITIONAL_GIS_ORDER)
获取投影​
# 获取栅格文件投影
ds: gdal.Dataset = gdal.Open(raster) # 打开栅格文件
prj: str = ds.GetProjection() # 获取投影信息
del ds # 释放内存
# 获取矢量文件投影
ds: ogr.DataSource = ogr.Open(vector) # 打开矢量文件
layer: ogr.Layer = ds.GetLayer() # 打开图层
spatialRef: str = layer.GetSpatialRef() # 获取投影信息
del ds, layer # 释放内存

CoordinateTransformation​

地理坐标转投影坐标​

from osgeo import osr
# 创建投影坐标系和地理坐标系
prosrs: osr.SpatialReference = osr.SpatialReference()
prosrs.ImportFromWkt(prj)
geosrs = prosrs.CloneGeogCS()
# 设置投影策略
prosrs.SetAxisMappingStrategy(osr.OAMS_TRADITIONAL_GIS_ORDER)
geosrs.SetAxisMappingStrategy(osr.OAMS_TRADITIONAL_GIS_ORDER)
# 建立投影变换矩阵
ct:osr.CoordinateTransformation = osr.CoordinateTransformation(geosrs, prosrs)
# 返回变换后坐标
coords: tuple = ct.TransformPoint(lon, lat)

投影坐标转地理坐标​

from osgeo import osr
# 读取或创建投影坐标系和地理坐标系
prosrs: osr.SpatialReference = osr.SpatialReference()
prosrs.ImportFromWkt(prj)
geosrs = prosrs.CloneGeogCS()
# 设置坐标策略
prosrs.SetAxisMappingStrategy(osr.OAMS_TRADITIONAL_GIS_ORDER)
geosrs.SetAxisMappingStrategy(osr.OAMS_TRADITIONAL_GIS_ORDER)
# 建立投影变换矩阵
ct:osr.CoordinateTransformation = osr.CoordinateTransformation(prosrs, geosrs)
# 返回变换后坐标
coords: tuple = ct.TransformPoint(x, y)

最佳实践​

矢量文件重投影​

from osgeo import osr
from osgeo import ogr
# 打开 shp
driver: ogr.Driver = ogr.GetDriverByName("ESRI Shapefile")
src_ds: ogr.DataSource = driver.Open("./cn_adm_pronvince/cn_adm_pronvince.shp")
# 创建投影转换
src_layer: ogr.Layer = src_ds.GetLayer()
src_prj: osr.SpatialReference = src_layer.GetSpatialRef()
dst_prj: osr.SpatialReference = osr.SpatialReference()
# 设置投影策略
src_prj.SetAxisMappingStrategy(osr.OAMS_TRADITIONAL_GIS_ORDER)
dst_prj.SetAxisMappingStrategy(osr.OAMS_TRADITIONAL_GIS_ORDER)
dst_prj.ImportFromEPSG(3857)
ct: osr.CoordinateTransformation = osr.CoordinateTransformation(
    src_prj, dst_prj)
# 创建 shp
dst_ds: ogr.DataSource = driver.CreateDataSource('./output/test.shp')
dst_layer: ogr.Layer = dst_ds.CreateLayer(
    'test', srs=dst_prj, geom_type=ogr.wkbPolygon, options=["ENCODING=UTF-8"])
# 创建字段
src_FeatureDefn: ogr.FeatureDefn = src_layer.GetLayerDefn()
for i in range(0, src_FeatureDefn.GetFieldCount()):
    dst_layer.CreateField(src_FeatureDefn.GetFieldDefn(i))
src_feature: ogr.Feature = src_layer.GetNextFeature()
dst_featureDefn: ogr.FeatureDefn = dst_layer.GetLayerDefn()
while (src_feature):
    dst_feature: ogr.Feature = ogr.Feature(dst_featureDefn)
    # 添加字段
    dst_fieldNum: int = dst_feature.GetFieldCount()
    for i in range(0, dst_fieldNum):
        dst_feature.SetField(i, src_feature.GetField(i))
    # 添加 geom
    geom: ogr.Geometry = src_feature.GetGeometryRef()
    geom.Transform(ct)
    dst_feature.SetGeometry(geom)
    # 添加 feature
    dst_layer.CreateFeature(dst_feature)
    src_feature = src_layer.GetNextFeature()
    del dst_feature

# 释放内存
del src_ds, dst_ds