python gdal.RasterizeLayer根据矢量创建栅格实现如下,RasterizeLayer支持shp、geojson等获取的矢量图层。
注意:burn_values=[1] 参数一定是序列[1],不能是1,否则报错:TypeError: not a sequence
import gdal
from osgeo import ogr
fg = 'polygon.geojson' # polygon.shp
fi = 'input.tif'
fl = 'output.tif'
# open origin image
data = gdal.Open(fi, gdal.GA_ReadOnly)
geo_transform = data.GetGeoTransform()
x_res = data.RasterXSize
y_res = data.RasterYSize
# set out iamge
y_ds = gdal.GetDriverByName('GTiff').Create(fl, x_res, y_res, 1, gdal.GDT_Byte,
options=['COMPRESS=LZW', 'BIGTIFF=YES'])
y_ds.SetGeoTransform(geo_transform)
# set out memory
target_ds = gdal.GetDriverByName('MEM').Create('', x_res, y_res, 1, gdal.GDT_Byte)
target_ds.SetGeoTransform(geo_transform)
band = target_ds.GetRasterBand(1)
band.SetNoDataValue(0)
# open vector
mb_v = ogr.Open(fg)
mb_l = mb_v.GetLayer()
if len(mb_l) > 0:
gdal.RasterizeLayer(target_ds, [1], mb_l, burn_values=[1])
y_buffer = band.ReadAsArray()
y_ds.WriteRaster(0, 0, x_res, y_res, y_buffer.tostring())
target_ds = None
y_ds = None
gdal