Rasterio: reproject not taking ndarray indexing into account

Created on 4 Feb 2018  路  5Comments  路  Source: mapbox/rasterio

Expected behavior and actual behavior.

In rasterio.warp.reproject, if the source ndarray (b) is indexing another ndarray (a), the reprojection should operate on b. It appears that is operates on a.

Steps to reproduce the problem.

a = np.zeros((5, 5))
b = a[::2, ::2]
reproject(
        b,
        destination,
        src_transform=src_transform,
        src_crs=src_crs,
        dst_transform=dst_transform,
        dst_crs=dst_crs,
        resampling=RESAMPLING.nearest)

Operating system

Linux version 4.10.0-42-generic (buildd@lgw01-amd64-012) (gcc version 6.3.0 20170406 (Ubuntu 6.3.0-12ubuntu2) )

Rasterio version and provenance

0.36.0 installed with conda install rasterio

bug

All 5 comments

So it looks like the ndarray raw data is passed to GDAL here:
https://github.com/mapbox/rasterio/blob/master/rasterio/_shim1.pyx#L109
But np.PyArray_DATA is the pointer to the base array data (a), so it is normal that b's indexing is not taken into account.

A simple fix would be to make a copy of the source ndarray, or at least to warn the user that array views are not taken into account.

I'm preparing a PR.

Thanks for looking into this @davidbrochart ! I think we likely have bugs like this in a couple other places, so a general solution would be very appreciated.

@davidbrochart I took another look at this issue and your PR this morning and realized that I need to ask you if switching to rasterio==1.0a12 addresses the problem. In recent versions we have switched to using GDAL's GDALRasterIOEx() function for its support of array (or view) striding: https://github.com/mapbox/rasterio/blob/master/rasterio/shim_rasterioex.pxi#L42-L43. I have a test of it at https://github.com/mapbox/rasterio/blob/master/tests/test_warp.py#L284.

Would you be willing to create a new conda environment and grab rasterio==1.0a12 from the conda-forge dev channel?

Yes, version 1.0a12 fixes the problem, no copy is needed anymore.
Thanks, this PR can be discarded.

Was this page helpful?
0 / 5 - 0 ratings

Related issues

ungarj picture ungarj  路  4Comments

snowman2 picture snowman2  路  3Comments

vincentsarago picture vincentsarago  路  4Comments

sgillies picture sgillies  路  5Comments

lwasser picture lwasser  路  3Comments