区域统计批处理功能（Python代码）

Tips ：本代码实现ArcGIS中的区域统计功能（Zonal），将统计的属性值直接赋予矢量文件属性表，可实现批量统计。

# -*- coding: utf-8 -*-
"""
Created on Wed Dec 16 16:13:59 2020

@author: xiao_gf

"""
import ogr, os
from rasterstats import zonal_stats

def zonal(input_vector,in_raster,Stats_type,raster_name):
            
        shp = ogr.Open(input_vector,1)
        lyr = shp.GetLayer() 
        # 添加字段
        zonal_Field = ogr.FieldDefn("{}".format(raster_name), ogr.OFTReal)
        # 添加字段名称为栅格文件名称，可根据自己需求修改
        zonal_Field.SetPrecision(9)
        lyr.CreateField(zonal_Field)
            
        zonal_method = zonal_stats(input_vector, in_raster, stats=[Stats_type])
        FID = 0
        for feat in lyr:
            Index = zonal_method[FID][Stats_type]
            feat.SetField("{}".format(raster_name), Index)
            feat.SetFID(FID)
            lyr.SetFeature(feat)
            FID +=1
        
if __name__ =="__main__":

    input_shp = r"D:\DATA\test.shp"
    in_folder = r"D:\DATA\raster"
    Stats_type = 'mean'
    # Stats_type表示统计的类型，包括：'min', 'max', 'mean', 'count', 'sum', 'std', 'median', 'majority', 'minority', 'unique', 'range
    files = os.listdir(in_folder)
    for file in files:
        if file[-4:] == '.tif':
            input_raster = os.path.join(in_folder,file)
            raster_name = file[:-4]
            zonal(input_shp, input_raster, Stats_type, raster_name)