使用rasterstats的zonal_stats遇width and height must be >0错误求助

我需要用Python的rasterstats包中的zonal_stats函数，批量从.tif栅格文件提取.shp矢量文件中每个面要素的栅格统计值（QGIS中可正常执行，但需处理200+文件，故需自动化）。

初始代码 #

import rasterio
import geopandas as gpd
import numpy as np
from rasterio.plot import show
from rasterstats import zonal_stats
from rasterio.transform import Affine


# 导入栅格文件
raster = rasterio.open(r'M:\PUBLIC\Felipe Dias\Pesquisa\Interpolação Espacial\Arroz_2019-03.tif')

# 读取栅格数值
array = raster.read(1)

# 获取仿射变换参数
affine = raster.transform

# 导入矢量文件
shapefile =  gpd.read_file(r'M:\PUBLIC\Felipe Dias\Pesquisa\Interpolação Espacial\Setores_Censit_SP_WGS84.shp')

# 计算分区统计
zs_shapefile = zonal_stats(shapefile, array, affine = affine,
                        stats=['min', 'max', 'mean', 'median', 'majority'])

运行报错 #

Input In [1] in <cell line: 22>
    zs_shapefile = zonal_stats(shapefile, array, affine = affine,

  File ~\Anaconda3\lib\site-packages\rasterstats\main.py:32 in zonal_stats
    return list(gen_zonal_stats(*args, **kwargs))

  File ~\Anaconda3\lib\site-packages\rasterstats\main.py:164 in gen_zonal_stats
    rv_array = rasterize_geom(geom, like=fsrc, all_touched=all_touched)

  File ~\Anaconda3\lib\site-packages\rasterstats\utils.py:41 in rasterize_geom
    rv_array = features.rasterize(

  File ~\Anaconda3\lib\site-packages\rasterio\env.py:387 in wrapper
    return f(*args, **kwds)

  File ~\Anaconda3\lib\site-packages\rasterio\features.py:353 in rasterize
    raise ValueError("width and height must be > 0")

尝试的解决方案及新问题 #

我尝试修改栅格的Affine符号，代码如下：

''' 尝试修改Affine参数 '''

old_tif = rasterio.open(r'M:\PUBLIC\Felipe Dias\Pesquisa\Interpolação Espacial\Arroz_2019-03.tif')
print(old_tif.profile) # 复制输出并修改符号

new_tif_profile = old_tif.profile

# 原Affine参数:
# Affine(0.004611149999999995, 0.0, -46.828504575,
#        0.0, 0.006521380000000008, -24.01169169)


new_tif_profile['transform'] = Affine(0.004611149999999995, 0.0, -46.828504575,
                                      0.0, -0.006521380000000008, 24.01169169)

new_tif_array = old_tif.read(1)
new_tif_array = np.fliplr(np.flip(new_tif_array))

with rasterio.open(r'M:\PUBLIC\Felipe Dias\Pesquisa\Interpolação Espacial\tentativa.tif', "w", **new_tif_profile) as dest:
    dest.write(new_tif_array, indexes=1)


dem = rasterio.open(r'M:\PUBLIC\Felipe Dias\Pesquisa\Interpolação Espacial\tentativa.tif')

# 读取栅格数值
array = dem.read(1)

# 获取仿射变换参数
affine = dem.transform

# 导入矢量文件
shapefile =  gpd.read_file(r'M:\PUBLIC\Felipe Dias\Pesquisa\Interpolação Espacial\Setores_Censit_SP_WGS84.shp')

# 计算分区统计
zs_shapefile = zonal_stats(shapefile, array, affine=affine,
                        stats=['min', 'max', 'mean', 'median', 'majority'])

修改后不再出现"width and height must be >0"错误，但zs_shapefile中的所有统计值均为None，无法满足需求。

请问该错误的原因是什么？我需要修改哪个符号才能解决问题？

内容的提问来源于stack exchange，提问作者femdias