Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

gdalwarp incorrectly re-projects data in EPSG:27700 to EPSG:3857 #3695

Closed
millionpoundhat opened this issue Apr 19, 2021 · 3 comments
Closed
Assignees

Comments

@millionpoundhat
Copy link

Expected behavior and actual behavior.

gdalwarp incorrectly re-projects data in EPSG:27700 to EPSG:3857.

There is 100m offset when re-projecting data from British National Grid to Web Mercator.

I believe the problem was introduced with Proj-8.

Steps to reproduce the problem.

Take the following WMS source xml as an example data source in 27700

  <Service name="TMS">
    <ServerUrl>https://api.os.uk/maps/raster/v1/wmts?key=yOMuONZcAGC54eFPLobKBAEmPGCIq3Gm&amp;service=WMTS&amp;request=GetTile&amp;version=2.0.0&amp;layer=Leisure_27700&amp;style=default&amp;outputFormat=image/png&amp;tileMatrixSet=EPSG:27700&amp;tileMatrix=9&amp;tileRow=${y}&amp;tileCol=${x}</ServerUrl>

  </Service>
  <DataWindow>
    <UpperLeftX>-238375</UpperLeftX>
    <UpperLeftY>1376256</UpperLeftY>
    <LowerRightX>700185.000,</LowerRightX>
    <LowerRightY>0</LowerRightY>
    <TileLevel>0</TileLevel>
    <TileX>0</TileX>
    <TileY>0</TileY>
    <SizeX>536320</SizeX>
    <SizeY>786432</SizeY>
    <YOrigin>top</YOrigin>
  </DataWindow>
  <Projection>EPSG:27700</Projection>
  <BlockSizeX>256</BlockSizeX>
  <BlockSizeY>256</BlockSizeY>
  <BandsCount>4</BandsCount>
  <MaxConnections>100</MaxConnections>
  <ZeroBlockHttpCodes>204,404</ZeroBlockHttpCodes>
</GDAL_WMS>

To test with

  • GDAL 3.2.2, released 2021/03/05
  • Proj Rel. 8.0.0, March 1st, 2021

I run the following

docker run -it -v "$(pwd)":/app osgeo/gdal:alpine-small-3.2.2 /bin/sh 

Then inside the container /app

gdalwarp -te 360443.864631 174074.392422 360828.386103 174457.931515 -te_srs epsg:27700 -t_srs epsg:3857 gdal-osgb-tms.xml gdal-3857-from-tms-3.2.2-gdal.tif

reprojects the tile incorrectly


To test with

  • GDAL 3.2.1, released 2020/12/29
  • Proj Rel. 7.2.1, January 1st, 2021

I ran the following

docker run -it -v "$(pwd)":/app osgeo/gdal:alpine-small-3.2.1 /bin/sh 

Then inside the container /app

gdalwarp -te 360443.864631 174074.392422 360828.386103 174457.931515 -te_srs epsg:27700 -t_srs epsg:3857 gdal-osgb-tms.xml gdal-3857-from-tms-3.2.1-gdal.tif

reprojects the tile correctly


I also tested on my Debian system, with the ubuntugis-stable PPA

  • Pop!_OS 20.10

  • Release: 20.10

  • Codename: groovy

  • GDAL 3.1.3, released 2020/09/01

  • Proj Rel. 7.1.0, August 1st, 2020

gdalwarp -te 360443.864631 174074.392422 360828.386103 174457.931515 -te_srs epsg:27700 -t_srs epsg:3857 gdal-osgb-tms.xml gdal-3857-from-tms-3.2.1-gdal.tif

I've also tested this with a geotiff in epsg:27700 incase there was an issue in the WMS/TMS driver and the error in reprojection remains.


Testing cs2cs

I also tested proj with cs2cs to see if the issue was there. Interestingly

Proj Rel. 7.1.0, August 1st, 2020
GDAL 3.1.3, released 2020/09/01

cs2cs EPSG:27700 EPSG:3857  <<EOF
360443.864631 174074.392422
EOF

Result = -286180.06 6703828.59 0.00

Proj Rel. 8.0.0, March 1st, 2021
GDAL 3.2.1, released 2020/12/29

cs2cs EPSG:27700 EPSG:3857  <<EOF
360443.864631 174074.392422
EOF

Result = -286180.20 6703828.38 0.00

There is a difference here but not 100m, I do note that the mercator projection was made more accurate in Proj 8 which could account for this.

@rouault rouault self-assigned this Apr 19, 2021
@rouault
Copy link
Member

rouault commented Apr 19, 2021

Confirmed

rouault added a commit to rouault/PROJ that referenced this issue Apr 19, 2021
…TERION_PARTIAL_INTERSECTION if area is specified"

This reverts commit ebe3425.

It was found to break gdalwarp usage in
OSGeo/gdal#3695 when passing a bbox that is
quite large.
@rouault
Copy link
Member

rouault commented Apr 19, 2021

This is indeed a PROJ 8.0.0 issue. Fix OSGeo/PROJ#2678 will be in PROJ 8.0.1

@rouault rouault closed this as completed Apr 19, 2021
@millionpoundhat
Copy link
Author

This is indeed a PROJ 8.0.0 issue. Fix OSGeo/PROJ#2678 will be in PROJ 8.0.1

Thanks!

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants