generating zonal stats from raster data
Note:This module is a thin layer on top of the rasterstats package to make its interface more compatible with the other geowrangler modules (e.g. vector zonal stats)
import matplotlib.pyplot as plt
import numpy as np
import rasterio
with rasterio.open(terrain_file) as src:
data = src.read(1)
data_crs = src.crs
data_bounds = src.bounds
print(data.shape, data_crs, data_bounds)
simple_aoi.total_bounds
ax = plt.imshow(data, cmap="pink")
ax = simple_aoi.plot(facecolor="none", edgecolor="blue")
%%time
results = create_raster_zonal_stats(
simple_aoi,
terrain_file,
aggregation=dict(func=["mean", "max", "min", "std"], column="elevation"),
extra_args=dict(nodata=np.nan),
)
results
Check that create_raster_zonal_stats uses the nodata attribute set in the tif file if extra_args.nodata is set to None
phl_adm = gpd.read_file("../data/region3_admin.geojson")
If the nodata parameter is explicitly set to -999999, the Population count should be > 0.
grid_aoi_results = create_raster_zonal_stats(
phl_adm,
"../data/phl_ppp_2020_constrained.tif",
aggregation=dict(
func=["sum"],
column="population",
output=["population_count"],
fillna=[True],
),
extra_args=dict(nodata=-99999),
)
grid_aoi_results
If the nodata parameter is not set, create_raster_zonal_stats should use the nodata attribute set by the geotiff file so that population count should still be > 0.
grid_aoi_results = create_raster_zonal_stats(
phl_adm,
"../data/phl_ppp_2020_constrained.tif",
aggregation=dict(
func=["sum"],
column="population",
output=["population_count"],
fillna=[True],
),
)
grid_aoi_results