Hey all,
I have some behaviour that I can’t wrap my head around. When I reproject a tif
using gdalwarp, the location (as visible in QGIS or after using gdal2tiles on a
map) changes. It is visible when checking the bounds using QGIS.
It is also noticeable when checking the bounds with rasterio:
```
orig.tif
EPSG:27573
bounds [828662.0710931281, 3203193.9700861075, 829595.8710931281,
3204344.6700861077]
EPSG:27573 bounds [828662.0710931281, 3203193.9700861075, 829595.8710931281,
3204344.6700861077]
EPSG:3857 bounds [578066.9358574335, 5479808.142856716, 579420.0775091586,
5481457.907968456]
EPSG:4326 bounds [5.192863637248721, 44.09275933293886, 5.205019115521942,
44.103402376252255] <----- 1
gdalwarp -t_srs EPSG:3857 orig.tif repro.tif
repro.tif
EPSG:3857
bounds [608015.5294828439, 6272484.913981743, 609507.0864301398,
6274297.75045133]
EPSG:27573 bounds [829484.0038619096, 3748101.0103735547, 830512.7970793804,
3749330.248761157]
EPSG:3857 bounds [608015.5294828439, 6272484.913981743, 609507.0864301398,
6274297.75045133]
EPSG:4326 bounds [5.46189643116462, 48.98599229672266, 5.475295315193526,
48.99667806803407] <----- 2
gdalwarp -t_srs EPSG:27573 repro.tif back.tif
back.tif
EPSG:27573
bounds [829484.0038619096, 3748101.034601565, 830512.7912119445,
3749330.2487611575]
EPSG:27573 bounds [829484.0038619096, 3748101.034601565, 830512.7912119445,
3749330.2487611575]
EPSG:3857 bounds [607947.003095218, 6272428.113526381, 609575.9130150724,
6274354.58877996]
EPSG:4326 bounds [5.461280848150926, 48.98565744919683, 5.475913594925499,
48.99701306471833] <----- 3
```
The accompanying python code is three times like the following
```
input_file = folder + 'orig.tif'
img = rasterio.open(input_file)
print('orig.tif')
print(img.crs)
print('bounds', list(warp.transform_bounds(img.crs, img.crs, *img.bounds)))
print('EPSG:27573 bounds',list(warp.transform_bounds(img.crs, 'EPSG:27573',
*img.bounds)))
print('EPSG:3857 bounds',list(warp.transform_bounds(img.crs, 'EPSG:3857',
*img.bounds)))
print('EPSG:4326 bounds', list(warp.transform_bounds(img.crs, 'EPSG:4326',
*img.bounds)))
```
As the arrows show, the bounds change between 1 and 2 (this is the unexpected
behaviour for me). Yet they don’t change when projecting back.
I am unsure if this is a bug in GDAL, something weird with this specific EPSG
or the source tif. Any further steps for investigation or ideas would be very
welcome
Many thanks
Evert
_______________________________________________
gdal-dev mailing list
[email protected]
https://lists.osgeo.org/mailman/listinfo/gdal-dev
