需求:有一个点shp文件和一个栅格,想要构建shp中每个点的缓冲区,并且缓冲区范围内的栅格值重新赋为0并输出新的tif文件
解决方法:使用python中的geopandas和rasterio中的掩膜操作实现
代码如下:
import numpy as np
import rasterio
from rasterio.features import geometry_mask
import geopandas as gpd
from shapely.geometry import Point
from shapely.ops import cascaded_union
from shapely.ops import unary_union
# 读取点数据
difangdef = gpd.read_file('abc.shp')
# 创建缓冲区
buffer_distance = 10000 # 缓冲区距离(单位:米)
buffers = difangdef.geometry.to_crs('EPSG:3395').buffer(buffer_distance)
# 将缓冲区合并为一个多边形
union_buffer = unary_union(buffers)
# 读取栅格数据
with rasterio.open('ddd.tif') as src:
# 利用栅格数据的空间参考创建一个蒙版
mask = geometry_mask([union_buffer], out_shape=src.shape, transform=src.transform, invert=True)
# 读取栅格数据作为numpy数组
raster_data = src.read(1)
# 将蒙版应用到栅格数据,这里的值可以任意修改赋的值,这里替换为-100
raster_data[mask] = -100
raster_data = np.where(raster_data == raster_data[0,0], np.nan, raster_data)
plt.imshow(raster_data)
plt.colorbar()
plt.show()
# 创建一个新的栅格数据文件来保存结果
with rasterio.open('path_to_output_raster.tif', 'w', **src.meta) as dst:
dst.write(raster_data, 1)
如果已有的是面文件,想要将面文件范围内的栅格值重新赋值,可以直接读取shp后进行掩膜操作,不需要构建缓冲区。文章来源:https://www.toymoban.com/news/detail-643002.html
上述代码实现后的结果如下:
文章来源地址https://www.toymoban.com/news/detail-643002.html
更多内容请关注公众号:GISerQ
到了这里,关于python:使用geopandas和rasterio将矢量范围内的栅格值赋为0并重新输出的文章就介绍完了。如果您还想了解更多内容,请在右上角搜索TOY模板网以前的文章或继续浏览下面的相关文章,希望大家以后多多支持TOY模板网!