gdal矢量转栅格_gdal 矢量数据转栅格

gdal矢量转栅格

def shp_to_tiff(shp_file, output_tiff, attribute):
    """

    :param shp_file:
    :param output_tiff:
    :param attribute: 定义栅格值的矢量属性
    :return:
    """
    start_time = datetime.datetime.now()
    print("start :" + str(start_time))
    # 读取shp文件
    driver = ogr.GetDriverByName("ESRI Shapefile")
    data_source = driver.Open(shp_file, 1)
    # 获取图层文件对象
    shp_layer = data_source.GetLayer()
    lon_min, lon_max, lat_min, lat_max = shp_layer.GetExtent()
    s_projection = str(shp_layer.GetSpatialRef())

    # (0,0,:,0,:,0)表示旋转系数
    # 自定义仿射矩阵系数 ， 1表示分辨率大小，决定了栅格像元的大小
    dst_transform = (lon_min, 1, 0, lat_max, 0, -1)
    d_lon = int(abs((lon_max - lon_min) / dst_transform[1])) # 除以横向分辨率
    d_lat = int(abs((lat_max - lat_min) / dst_transform[5]))

    # 根据模板tif属性信息创建对应标准的目标栅格
    target_ds = gdal.GetDriverByName('GTiff').Create(output_tiff, d_lon, d_lat, 1, gdal.GDT_Byte)
    target_ds.SetGeoTransform(dst_transform)
    target_ds.SetProjection(s_projection)

    band = target_ds.GetRasterBand(1)
    # 设置背景数值
    NoData_value = 0
    band.SetNoDataValue(NoData_value)
    band.FlushCache()

    # 调用栅格化函数。gdal.RasterizeLayer函数有四个参数，分别有栅格对象，波段，矢量对象，value的属性值将为栅格值
    option = 'ATTRIBUTE=%d' % (attribute)
    gdal.RasterizeLayer(target_ds, [1], shp_layer, options=[option])
    # 直接写入
    y_buffer = band.ReadAsArray()
    target_ds.WriteRaster(0, 0, d_lon, d_lat, y_buffer.tobytes())
    start_time = datetime.datetime.now()
    print("end :" + str(start_time))
    target_ds = None  # todo 释放内存，只有强制为None才可以释放干净
    del target_ds, shp_layer
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46