基于栅格单元与栅格单元的分区统计

 输入的区域可以是栅格的或者是矢量的，值数据只能是栅格的

import numpy as np
import pandas as pd


def zonal_statistic(zonal_arr, data_arr, stats_field=['mean'], nodata=255):
    """
    Parameters
    ----------
    zonal_arr:分区统计的区域数据
    data_arr:分区统计中需要作统计的数据
    stats_field:分区统计的字段，有均值、最大最小值、总数等
    nodata: 区域中nodata值，默认为255
    :rtype: object
    """
    classes_list = np.array(sorted(np.unique(zonal_arr)))
    classes_list = classes_list[classes_list != nodata]
    classes_list = classes_list[~np.isnan(classes_list)]
    if classes_list.size == 0:
        return pd.DataFrame()
    zonal_sta_list = []
    zonal_sta = pd.Series()
    for i in classes_list:
        if 'mean' in stats_field:
            zonal_sta['mean'] = np.nanmean(data_arr[zonal_arr == i])
        else:
            zonal_sta['mean'] = np.nan
        if 'max' in stats_field:
            zonal_sta['max'] = np.nanmax(data_arr[zonal_arr == i])
        else:
            zonal_sta['max'] = np.nan
        if 'min' in stats_field:
            zonal_sta['min'] = np.nanmin(data_arr[zonal_arr == i])
        else:
            zonal_sta['min'] = np.nan
        if 'count' in stats_field:
            zonal_sta['count'] = len(data_arr[zonal_arr == i])
        else:
            zonal_sta['count'] = np.nan
        zonal_sta_copy = zonal_sta.copy()
        zonal_sta_list.append(zonal_sta_copy)

    return pd.DataFrame(zonal_sta_list, index=classes_list)[stats_field]