[GIS] gdal2tiles misalignment

gdalgdal2tiles

I have a map image that I have reprojected from OSGB36 to googleMaps projection (EPSG:900913). When I run gdal2tiles on the reprojected map, the results are not alligned with the google maps images (its not just out by a little bit but by 10's of miles).

If I run it through MapTiler (from maptiler.org) it is aligned fine. MapTiler asks for an SRS which I dont specify when using gdal2tiles. So I'm thinking I've missied out a stage?

Also when I run gdalinfo on the warped file, the coords of the corners seem wrong. E.g. for the bottom-left corner the latitude is 51.143906, but it should be 50.955887 for google maps.

I'm usign gdal16 installed by OSGeo4W on Windows7 64-bit.

Warped file here..
https://docs.google.com/leaf?id=0B93lksnTC7_cNzIzZWQ0MDYtMTY5Yy00NjY0LThiMjktZjc2YzA3MmU5NjI4

gdalinfo on warped file…

C:\OSGeo4W\apps\gdal-16\bin>gdalinfo e:buckland_WARPED.tiff
Driver: GTiff/GeoTIFF
Files: e:buckland_WARPED.tiff
Size is 611, 419
Coordinate System is:
PROJCS["Google Maps Global Mercator",
    GEOGCS["WGS 84",
        DATUM["WGS_1984",
            SPHEROID["WGS 84",6378137,298.2572235630016,
                AUTHORITY["EPSG","7030"]],
            AUTHORITY["EPSG","6326"]],
        PRIMEM["Greenwich",0],
        UNIT["degree",0.0174532925199433],
        AUTHORITY["EPSG","4326"]],
    PROJECTION["Mercator_1SP"],
    PARAMETER["central_meridian",0],
    PARAMETER["scale_factor",1],
    PARAMETER["false_easting",0],
    PARAMETER["false_northing",0],
    UNIT["metre",1,
        AUTHORITY["EPSG","9001"]]]
Origin = (-474051.292871501410000,6617009.987366309400000)
Pixel Size = (8.392224865731517,-8.392224865731517)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left  ( -474051.293, 6617009.987) (  4d15'30.51"W, 51d 9'49.57"N)
Lower Left  ( -474051.293, 6613493.645) (  4d15'30.51"W, 51d 8'38.06"N)
Upper Right ( -468923.643, 6617009.987) (  4d12'44.69"W, 51d 9'49.57"N)
Lower Right ( -468923.643, 6613493.645) (  4d12'44.69"W, 51d 8'38.06"N)
Center      ( -471487.468, 6615251.816) (  4d14'7.60"W, 51d 9'13.82"N)
Band 1 Block=611x3 Type=Byte, ColorInterp=Red
  Mask Flags: PER_DATASET ALPHA
Band 2 Block=611x3 Type=Byte, ColorInterp=Green
  Mask Flags: PER_DATASET ALPHA
Band 3 Block=611x3 Type=Byte, ColorInterp=Blue
  Mask Flags: PER_DATASET ALPHA
Band 4 Block=611x3 Type=Byte, ColorInterp=Alpha

The sequence of commands I run are…

Georeference the OSGB PNG images…

gdal_translate -a_srs EPSG:27700 
-gcp 0 400 241538.5 119812 
-gcp 600 0 244713.5 121928.66666667 
-gcp 0 0 241538.5 121928.66666667 
-gcp 600 400 244713.5 11981
bucklands1.png bucklands1.tiff

Warp into googleMaps projection (EPSG:900913) …

gdalwarp -s_srs EPSG:27700 -t_srs EPSG:900913 -of GTiff
bucklands1.tiff bucklands1_warped.tiff

Tile (make into googleMap 256×256 tiles) …

gdal2tiles bucklands1_warped.tiff gmTiles

It works if I by-pass the gdalwarp stage and make google tiles directly from the OSGB36 images. But the image quality is poor as the interpolation is not so good.

Best Answer

It turns out that if you replace EPSG:900913 with either EPSG:3857 or EPSG:3785 in the gdal warping you get the google-map tiles correctly aligned. Apparently these codes are alternatives for the unofficial googleMaps EPSG:900913. Though there seems to be lots of confusion about which is the correct code - e.g. MapTiler doesnt recognise EPSG:3857.

Anyway this works...

gdalwarp -s_srs EPSG:27700 -t_srs EPSG:3857 -of GTiff
bucklands1.tiff bucklands1_warped.tiff
Related Question