前言

目前有一张tif格式的栅格影像,需要在web地图上进行展示,使用动态切片wms的方式,渲染速度比较慢,而且大的时候会出现模糊的问题。并且后面需要做多期影像的切换,渲染与加载效率也值得关注。

计划是使用栅格转矢量的方式,将栅格数据转为矢量shp文件,然后进行矢量切片,使用mapbox进行前端动态渲染。在网上查询了很多资料,有人说使用d3-contour在node.js中生成或者使用rasterio在python中进行转换,整体过程都比较麻烦,很不易实现。最终选定了使用gdal进行栅格转矢量的方法,代码比较简单。
原始tif影像(12.8mb)如下:

核心函数

gdal中栅格转矢量的函数主要是以下两个,二者的参数没有任何区别,只是功能有区别:

fpolygonize(*args, **kwargs)

fpolygonize(band srcband, band maskband, layer outlayer, int ipixvalfield, char options=none, gdalprogressfunc callback=0, void * callback_data=none) -> int

将每个像元转成一个矩形。

polygonize(*args, **kwargs) **

polygonize(band srcband, band maskband, layer outlayer, int ipixvalfield, char ** options=none, gdalprogressfunc callback=0, void * callback_data=none) -> int

将每个像元转成一个矩形,然后将相似的像元进行合并。

转换代码

from osgeo import gdal, ogr, osr
import os
import datetime
import numpy as np

path = "z_nafp20210727.tif"


if __name__ == '__main__':
    start_time = datetime.datetime.now()

    inraster = gdal.open(path)  # 读取路径中的栅格数据
    inband = inraster.getrasterband(1)  # 这个波段就是最后想要转为矢量的波段,如果是单波段数据的话那就都是1
    prj = osr.spatialreference()
    prj.importfromwkt(inraster.getprojection())  # 读取栅格数据的投影信息,用来为后面生成的矢量做准备

    outshp = path[:-4] + ".shp"  # 给后面生成的矢量准备一个输出文件名,这里就是把原栅格的文件名后缀名改成shp了
    drv = ogr.getdriverbyname("esri shapefile")
    if os.path.exists(outshp):  # 若文件已经存在,则删除它继续重新做一遍
        drv.deletedatasource(outshp)
    polygon = drv.createdatasource(outshp)  # 创建一个目标文件
    poly_layer = polygon.createlayer(path[:-4], srs=prj, geom_type=ogr.wkbmultipolygon)  # 对shp文件创建一个图层,定义为多个面类
    newfield = ogr.fielddefn('value', ogr.oftreal)  # 给目标shp文件添加一个字段,用来存储原始栅格的pixel value,浮点型,
    poly_layer.createfield(newfield)

    gdal.polygonize(inband, none, poly_layer, 0)  # 核心函数,执行的就是栅格转矢量操作
    # gdal.fpolygonize(inband, none, poly_layer, 0)  # 只转矩形,不合并
    polygon.synctodisk()
    polygon = none
    end_time = datetime.datetime.now()
    print("succeeded at", end_time)
    print("elapsed time:", end_time - start_time)  # 输出程序运行所需时间

转换效果

  • 使用fpolygonize

转换之后的矢量数据有270mb,非常大,打开非常卡

  • 使用polygonize

合并之后的矢量数据有48mb,相对第一种方法数据量大大减少

到此这篇关于python 使用gdal实现栅格tif转矢量shp的文章就介绍到这了,更多相关python栅格tif转矢量shp内容请搜索www.887551.com以前的文章或继续浏览下面的相关文章希望大家以后多多支持www.887551.com!