I wrote a Python script who calculate, for rasters, mean value for each polygon of a shapefile. But I think that it's a "war machine" (it takes many treatment time) and I would like to simplify it.

Is somebody know a simple way to do it?

I'm using geopandas for shapefile manipulations.


You can use the rasterstats package for this.

For example (assuming you have a geopandas.GeoDataFrame called gdf):

from rasterstats import zonal_stats

with rasterio.open("/path/to/raster.tif") as src:
    affine = src.transform
    array = src.read(1)
    df_zonal_stats = pd.DataFrame(zonal_stats(gdf, array, affine=affine))

# adding statistics back to original GeoDataFrame
gdf2 = pd.concat([gdf, df_zonal_stats], axis=1) 
  • It works perfectly! but I forgot that I want always more :D and I would like to have 2 more value: 5 and 95 percentiles... apparently it isn't possible. – Tim C. Sep 26 '18 at 15:12
  • 2
    You can pass custom defined functions, see pythonhosted.org/rasterstats/… – joris Sep 26 '18 at 19:07
  • Awesome, all I need ! – Tim C. Sep 27 '18 at 6:58

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.