基于python对栅格数据批量分区统计

基于python对栅格数据批量分区统计

前言

如何基于python对tif格式的栅格数据批量进行分区统计。使用 Python 对 .tif 格式的栅格数据进行批量分区统计，通常可以通过 rasterio 和 numpy 等库来实现。具体步骤包括读取栅格数据、进行分区处理（如基于掩膜或特定区域划分），然后计算统计信息（如均值、最大值、最小值等）。以下是实现这一流程的典型代码框架：

---

一、主要步骤

读取栅格文件：使用 rasterio 来加载 .tif 栅格数据。
定义分区掩膜：可以根据已有的矢量边界文件（如 .shp 文件）或自定义掩膜对栅格数据进行区域划分。
分区统计：对每个区域进行统计分析，如均值、最大值、最小值等。

二、使用步骤

1.代码

# -*- encoding:utf-8 -*-

import geopandas as gpd
import rasterio
from rasterio.mask import mask
import numpy as np
import pandas as pd
import os

# 读取shapefile文件
shp_path = r'shp数据.shp'
gdf = gpd.read_file(shp_path)

# 获取FID字段
fid_field = 'Ecosystem'  # 根据你的shapefile的FID字段名称进行修改

# 存放结果的列表
results = []

# 读取每年的tif文件夹路径
tif_folder = r'tif数据'
tif_files = [os.path.join(tif_folder, f) for f in os.listdir(tif_folder) if f.endswith('.tif')]

# 遍历所有tif文件
for tif_file in tif_files:
    year = os.path.basename(tif_file).split('.')[0]  # 假设文件名包含年份信息

    # 打开tif文件
    with rasterio.open(tif_file) as src:
        # 遍历每个FID对应的geometry
        for _, row in gdf.iterrows():
            geom = [row['geometry']]
            fid = row[fid_field]

            # 裁剪栅格数据
            out_image, out_transform = mask(src, geom, crop=True)
            out_image = out_image[0]  # 提取第一个波段

            # 计算非空像元的均值
            if np.any(out_image):
                mean_value = np.mean(out_image[out_image != src.nodata])
            else:
                mean_value = np.nan
#def batch_process_rasters(raster_dir, shp_path):
  #  results = []
   # for filename in os.listdir(raster_dir):
   #     if filename.endswith(".tif"):
      #      raster_path = os.path.join(raster_dir, filename)
       #     masked_raster = mask_raster_with_shapefile(raster_path, shp_path)
       #     stats = calculate_statistics(masked_raster)
         #   results.append({
         #       'file': filename,
          #      'mean': stats[0],
          #      'max': stats[1],
           #     'min': stats[2]
       #     })
   # return results


            # 保存结果
            results.append({
                'FID': fid,
                'Year': year,
                'Mean': mean_value
            })

# 将结果转换为DataFrame并保存为CSV
df = pd.DataFrame(results)
df.to_csv(r'F:\output.csv', index=False)

print("统计完成，结果已保存到output.csv")

---

解释：

read_raster 函数：用于读取 .tif 栅格文件的第一波段数据。
read_shapefile 函数：读取分区掩膜的矢量文件（例如 .shp 文件），用于后续的区域裁剪。
mask_raster_with_shapefile 函数：将矢量文件作为掩膜，裁剪栅格文件，从而对特定区域进行分析。
calculate_statistics 函数：针对每个裁剪后的栅格区域，计算统计信息（如均值、最大值、最小值）。
batch_process_rasters 函数：对多个 .tif 文件批量进行处理，并计算每个文件的分区统计值。

依赖库：

rasterio：用于读取和操作栅格数据。
geopandas：用于处理矢量数据（如 .shp 文件）。
numpy：用于计算统计信息。

注意事项：

NoData 值处理：代码中假设栅格文件存在 NoData 值（如 -9999），统计计算时需排除这些值。
分辨率匹配：确保栅格和矢量文件的坐标参考系相同，避免出现坐标系统不匹配的情况。
批量处理：通过遍历目录实现对多个 .tif 文件的批量处理。