0

I'm trying to overlay tiles from a WMS server (a map of the ocean biomass) on a Google Map (android), but the projection is slightly distorted.

As you can see in this picture, the longitudes are OK, but the latitudes get distorted the farther away from the equator you get:

Map distortion

When zooming in, the distortion is no longer visible:

enter image description here

The WMS accepts EPSG:4326 bounding boxes (I've reviewed other questions on this topic, but they all use a WMS that accepts different coordinate systems).

I'm calculating the WMS BBOX coordinates out of the (x,y,zoom) parameters in getTile by calling tile2boundingBox(x,y,zoom) (calculation taken from https://wiki.openstreetmap.org/wiki/Slippy_map_tilenames):

double tile2lon(int x, int z) {
  return x / pow(2.0, z) * 360.0 - 180;
}

double tile2lat(int y, int z) {
  double n = pi - (2.0 * pi * y) / pow(2.0, z);
  return vector_math.degrees(atan(sinh(n)));
}

tile2boundingBox(final int x, final int y, final int zoom) async {
  final north = tile2lat(y, zoom);
  final south = tile2lat(y + 1, zoom);
  final west = tile2lon(x, zoom);
  final east = tile2lon(x + 1, zoom);

  return [south, west, north, east];
}

And then I create a WMS url to fetch the tile with the BBOX e.g:

https://my.cmems-du.eu/thredds/wms/global-reanalysis-bio-001-029-daily?SERVICE=WMS&VERSION=1.3.0&REQUEST=GetMap&FORMAT=image%2Fpng&TRANSPARENT=true&LAYERS=chl&STYLES=boxfill%2Frainbow&TILED=true&COLORSCALERANGE=0.049000002%2C10.9050001&WIDTH=512&HEIGHT=512&CRS=EPSG%3A4326&BBOX=0.0,-90.0,40.97989806962013,-45.0

I've also compared my calculation to the python program in MapTiler: https://www.maptiler.com/google-maps-coordinates-tile-bounds-projection/ and the coordinated are exactly the same:

enter image description here

What am I doing wrong? Thanks.

3
  • 2
    Since Google Maps use projected EPSG:3857, you should use the same CRS also for your WMS service calls and you should be OK.
    – TomazicM
    Aug 10 '20 at 16:55
  • Looking at the GetCapabilities document from the linked service, it supports EPSG:3857, I tried pulling it up in QGIS, and overlaying it over an OSM layer both in EPSG:3857 and it appears to align.
    – Andy C
    Aug 10 '20 at 17:53
  • OK, so what's wrong with the calculations I pasted? What is the correct calculation to get from google maps tile (x,y, zoom) to the WMS BBOX?
    – syonip
    Aug 11 '20 at 11:11
1
+100

Firstly, I recommend using an existing library rather than implementing your own XYZ tile to latitude/longitude conversion. https://github.com/mapbox/mercantile is good if you're using Python. There's much less likely to be bugs in popular open source projects than in a roll-your-own.

To use this, you would:

import mercantile
mercantile.xy_bounds(2, 2, 3)
Bbox(left=-10018754.171394622, bottom=5009377.085697311, right=-5009377.085697311, top=10018754.171394622)

Substitute those into your WMS call, using CRS=EPSG:3857 as well:

https://my.cmems-du.eu/thredds/wms/global-reanalysis-bio-001-029-daily?SERVICE=WMS&VERSION=1.3.0&REQUEST=GetMap&FORMAT=image%2Fpng&TRANSPARENT=true&LAYERS=chl&STYLES=boxfill%2Frainbow&TILED=true&COLORSCALERANGE=0.049000002%2C10.9050001&WIDTH=512&HEIGHT=512&BBOX=-10018754.171394622,5009377.085697311,-5009377.085697311,10018754.171394622&CRS=EPSG:3857


Secondly, you shouldn't need to perform this conversion at all! What the comments on your post suggest is that the service supports using EPSG:3857 (Web Mercator) which is the typical projection of XYZ tile services.

Compare:

  1. https://my.cmems-du.eu/thredds/wms/global-reanalysis-bio-001-029-daily?SERVICE=WMS&VERSION=1.3.0&REQUEST=GetMap&FORMAT=image%2Fpng&TRANSPARENT=true&LAYERS=chl&STYLES=boxfill%2Frainbow&TILED=true&COLORSCALERANGE=0.049000002%2C10.9050001&WIDTH=512&HEIGHT=512&BBOX=-10948228.4353,2993885.5239,-4872401.9310,9147983.5452&CRS=EPSG:3857
  2. https://my.cmems-du.eu/thredds/wms/global-reanalysis-bio-001-029-daily?SERVICE=WMS&VERSION=1.3.0&REQUEST=GetMap&FORMAT=image%2Fpng&TRANSPARENT=true&LAYERS=chl&STYLES=boxfill%2Frainbow&TILED=true&COLORSCALERANGE=0.049000002%2C10.9050001&WIDTH=512&HEIGHT=512&BBOX=0.0,-90.0,40.0,-45.0&CRS=EPSG:4326

The former is Web Mercator (EPSG:3857), and the latter is WGS84 (EPSG:4326). However note that your BBOX parameter must match the coordinate system that you request. You cannot request an output in Web Mercator if you give it a WGS84 bounding box. (Well, you can try, but you won't get sensible output.)

3
  • Thanks, but as I said I've tested the calculation using a third party library map tiler, and I'm getting the same results. The second link you've attached doesn't work by the way. Let's try this: what is the correct WMS call for the tile: (x=2,y=2,z=3)? I'll manually overlay the image on top of the google map tile in photoshop and check if it fits.
    – syonip
    Aug 17 '20 at 16:06
  • @syonip see edits Aug 17 '20 at 21:19
  • Great thanks! I used a projection from EPSG:4326 to EPSG:3857 and it fits correctly. I still don't understand why sending the BBOX in WGS84 gives a distorted result, as the WMS seems to support both.
    – syonip
    Aug 18 '20 at 7:53

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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