5

I am working with the Rasterstats package to run zonal statistics on two raster datasets. The documentation calls for the following:

import rasterio
from rasterstats import zonal_stats

with rasterio.open('tests/data/slope.tif') as src:
    affine = src.affine
    array = src.read(1)
zs = zonal_stats('tests/data/polygons.shp', array, affine=affine)

However, src.affine is no longer supported in Rasterio. Instead, the documentation recommends using src.transform, which yields something like:

Affine(1.0, 0.0, 670625.0016573452,
       0.0, -1.0, 5282947.795216726)

When I do something like the following, I get an error:

from rasterstats import zonal_stats

dataset = rasterio.open('/path/to/image.tif')
affine = dataset.transform * (0,0) # upper left corner
zs = zonal_stats(segments, ndvi, affine=affine, stats="mean")

IndexError: only integers, slices (`:`), ellipsis (`...`), numpy.newaxis (`None`) and integer or boolean arrays are valid indices

What affine format does the Rasterstats zonal_stats function require?

0

1 Answer 1

4

The results of dataset.transform and dataset.affine are equals (see Migrating to Rasterio 1.0: affine.Affine() vs. GDAL-style geotransforms)

With rasterio version < 1 (affine = affine transform, transform = GDAL-style geotransforms,

import rasterio 
rasterio.__version__
'0.36.0'
with rasterio.open('slope.tif') as src:
   affine = src.affine
   transform = src.transform
affine
Affine(25.52514657450613, 0.0, 244300.61494985913,
   0.0, -25.52514657450613, 1000868.7876863468)
print(affine)
| 25.53, 0.00, 244300.61|
| 0.00,-25.53, 1000868.79|
| 0.00, 0.00, 1.00|
transform
[244300.61494985913, 25.52514657450613, 0.0, 1000868.7876863468, 0.0, -25.52514657450613]

And the conversion

from affine import Affine
# or
from rasterio import Affine
aff = Affine.from_gdal(*transform)
aff
Affine(25.52514657450613, 0.0, 244300.61494985913,
   0.0, -25.52514657450613, 1000868.7876863468)
gda = aff.to_gdal()
gda
(244300.61494985913, 25.52514657450613, 0.0, 1000868.7876863468, 0.0, -25.52514657450613)

With rasterio version > 1.0b (transform = affine transform)

import rasterio
rasterio.__version__
'1.0.13'
with rasterio.open('slope.tif') as src:
    transform = src.transform
transform
Affine(25.52514657450613, 0.0, 244300.61494985913,
   0.0, -25.52514657450613, 1000868.7876863468)
print(transform)
| 25.53, 0.00, 244300.61|
| 0.00,-25.53, 1000868.79|
| 0.00, 0.00, 1.00|
transform.to_gdal()
(244300.61494985913, 25.52514657450613, 0.0, 1000868.7876863468, 0.0, -25.52514657450613)

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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