python 使用GDAL实现栅格tif转矢量shp的方式小结
作者:GIS开发者 发布时间:2021-10-02 07:13:56
前言
目前有一张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,相对第一种方法数据量大大减少
来源:https://blog.csdn.net/GISuuser/article/details/119452042
猜你喜欢
- 以下为引用的内容: <html> <head> <title>不刷新页面查询的方法&
- 导言在前面的教程里我们学习了DataList提供了一些风格样式的属性.而且我们还学习了如何定义HeadStyle, ItemStyle, A
- 最近老婆大人的公司给老婆大人安排了一个根据关键词查询google网站排名的差事。老婆大人的公司是做seo的,查询的关键词及网站特别的多,看着
- RegExp对象的语义和使用:检查字符串匹配获取字符串中的部分内容在原字符串的基础上构建一个新的字符串(包括添加、删除和修改)构建一个Reg
- 最近在学习PHP,以下是看PHP100视频教程,做的学习笔记,在这里存放以便今后使用。apache--PHP--DB(mysql)一、apa
- 前言:通常我们创建类都是使用class 类名,但是小伙伴们有没有想过,类是由谁来创建的呢,python中常说的万物皆对象,对象是由类创建的,
- 这几天在落伍上转转,发现有朋友不太明白一些网站在会员注册时,当输入用户名后没按“确定”提交数据,系统也能马上检测该用户名是否已经存在。在此我
- 制作自己的训练集下图是我们数据的存放格式,在data目录下有验证集与测试集分别对应iris_test, iris_train 为了向伟大的M
- 代码如下: <% dim fso,objFolder,objFiles dim filelist Set fso=Server.Cre
- 在程序的开发过程中,处理分页是大家接触比较频繁的事件,因为现在软件基本上都是与数据库进行挂钓的。但效率又是我们所追求的,如果是像原来那样把所
- Event Handler在HDA中,要创建Python脚本,需要先选择一个事件处理器(EventHandle),他表示你要在什么时候执行你
- 在做我的友情链接批量检查工具过程中,碰到一些情况,就是对方网页会用gzip压缩。用gzip压缩的好处是,能压缩网页大小,加快网页的浏览速度,
- ASP具备动态输出任一Office应用程序文件格式的功能。在开始编写代码之前,我们首先需要做的就是设置正确的文件类型,因为浏览器需要知道如何
- 首先介绍下比较简单但必不可少且实用的知识,可以当手册查询,适合像我一样的新手看。PHP常用库函数介绍一、PHP字符串操作常用函数1.确定字符
- 一、使用场景需求1、在实际项目开发过程中,用户可以操作的数据,我们往往会新增一个字段,来保存用户最后一次修改时间2、一些系统中,我们需要存储
- SQL语言共分为四大类:数据查询语言DQL,数据操纵语言DML, 数据定义语言DDL,数据控制语言DCL。其中用于定义数据的结构,比如 创建
- Python脚本编辑使用Python对fasta格式的序列进行基本信息统计预期设计输出文件中包括fasta文件名,序列长度,GC含量以及AT
- Abs (数值)绝对值。一个数字的绝对值是它的正值。空字符串 (null) 的绝对值,也是空字符串。未初始化的变数,其绝对为 0例子:ABS
- 破解百度翻译翻译是一件麻烦的事情,如果可以写一个爬虫程序直接爬取百度翻译的翻译结果就好了,可当我打开百度翻译的页面,输入要翻译的词时突然发现
- 网上有许关于固定表格的标题行的文章,但是既要固定标题行又要固定标题列的却几乎没有。现我写下如下代码以供大家参考:<html> &