【Python&GIS】面矢量数据投影转换（WGS84转地方坐标系）

之前分享过点矢量怎么进行投影转换，今天跟大家分享下面矢量如何投影转换。代码与之前的类似，只要注意一下GDAL对矢量文件划分的关系层次即可。

ogr库是一个处理地理空间矢量数据的开源库。它可以读取多种数据格式，进行地理处理、属性表操作、数据分析等操作。目前ogr和osr库已集成到GDAL库中，可以对栅格数据、矢量数据进行处理分析，被3S的研究人员广泛应用。感兴趣的可以自己去了解一下，不懂得可以一起交流！

1.加载GDAL库（水字数）

from osgeo import ogr, osr

2.定义目标投影坐标系

target_proj = osr.SpatialReference()
# 初始化osr.SpatialReference对象以形成一个合法的坐标系统
target_proj.ImportFromEPSG(4326)
# 向对象中写入WGS84坐标系统

3.读取矢量中每个特征的投影信息

.GetSpatialRef()函数可以获取图层中的坐标系信息，注意该函数的作用空间为图层。图层用GetSpatialRef()获取坐标信息，元素用GetSpatialReference()获取坐标信息。

ds = ogr.Open(path)
layer = ds.GetLayer(0)
feature_count = layer.GetFeatureCount()
# 查看矢量包含多少个特征
for i in range(feature_count):
    # 循环所有特征
    """spatialRef_layer = layer.GetSpatialRef()
    print(spatialRef_layer)
    # 图层用GetSpatialRef()获取坐标信息，元素用GetSpatialReference()获取坐标信息"""
    feature = layer.GetFeature(i)
    geom = feature.GetGeometryRef()
    source_proj = geom.GetSpatialReference()
    # 获取特征的原始坐标系

4.坐标转换参数

第一个参数是数据的坐标系，第二个参数是目标坐标系。

transform = osr.CoordinateTransformation(source_proj, target_proj)
# 定义坐标系转换参数    
geom.Transform(transform)

5.完整代码

# -*- coding: utf-8 -*-
"""
@Time ： 2023/5/18 16:26
@Auth ： RS迷途小书童
@File ：Face Vector Data Projection Transform.py
@IDE ：PyCharm
@Purpose ：对面矢量文件进行投影转换
"""
from osgeo import ogr, osr
"""layer图层才能调用投影，数据资源是不能调用投影的。矢量数据分为datasource，layer，feature三个层次"""


def Projection_Transform(path):
    """
    :param path: 输入需要投影转换的面矢量
    :return:
    """
    target_proj = osr.SpatialReference()
    # 初始化osr.SpatialReference对象以形成一个合法的坐标系统
    target_proj.ImportFromEPSG(4326)
    # 向对象中写入WGS84坐标系统
    ds = ogr.Open(path)
    layer = ds.GetLayer(0)
    feature_count = layer.GetFeatureCount()
    # 查看矢量包含多少个特征
    for i in range(feature_count):
        # 循环所有特征
        """spatialRef_layer = layer.GetSpatialRef()
        print(spatialRef_layer)
        # 图层用GetSpatialRef()获取坐标信息，元素用GetSpatialReference()获取坐标信息"""
        feature = layer.GetFeature(i)
        geom = feature.GetGeometryRef()
        source_proj = geom.GetSpatialReference()
        # 获取特征的原始坐标系
        transform = osr.CoordinateTransformation(source_proj, target_proj)
        # 定义坐标系转换参数
        geom.Transform(transform)
        # print(geom.GetSpatialReference())
        ds = None
    del ds


if __name__ == "__main__":
    input_path = 'B:/1.shp'
    Projection_Transform(input_path)