Google Web Mercator: ambiguous coordinate system

On October 1, 2014, the US National Geospatial Intelligence Agency ( NGA ) published a report outlining the criticism of the Web Mercator coordinate system used in many web-based map services. A detailed explanation of the problem and recommendations for NGA partners were attached to the document. The document received a great response, but not all articles based on this report were distinguished by their accuracy and literacy. This concerns, for example, an article on the GIS Association website, which, due to gross errors in terminology, can be considered illiterate. Since it is with this coordinate system that web services developers encounter most often, I think it makes sense to sort out the problem.

What is it about?

To begin with - a couple of definitions, without which some details can not be clear. It is important to understand that Web Mercator is a coordinate system, and not just a projection, although its name resembles the well-known Mercator projection. This terminological discrepancy misleads the readers of the article on the GIS Association website. The difference between the projection and the coordinate system is that the projection is only the way in which the complex shape of the Earth's figure model turns into a plane, while the coordinate system also includes the mathematical definition of a model (ellipsoid or spheroid) that approximates the complex shape of the Earth.
')


In this illustration, red indicates that it refers only to the projection mechanism (in this case, cylindrical). The coordinate system also refers to everything that is shown here.

In turn, it is this very approximating model of the surface (the dotted sphere in the figure above, which defines the coordinates λ, φ ) and is the source of the problem, which will be discussed further.

Story

I can not say for certain who the first and when it all came to mind. But, as far as I know, the first major project that began using the Web Mercator coordinate system was Google Maps service, and it happened in 2005. The developers then had the task to simplify the calculations necessary for working with cartographic data, and the most obvious thing that could be done was to use a sphere in the coordinate system instead of an ellipsoid. It is interesting that Gerard Mercator himself most likely proceeded from the same geometric representations, creating his own way of projecting maps onto a plane, because only Newton, who lived a little later, suggested the hypothesis that the Earth, due to centrifugal force, has the shape of an ellipsoid of rotation and not the ball. Thus, the developers of Google, in a sense, returned to the sixteenth century.

Criticism of this approach in professional circles does not sound the first time. Since 2005, the European Petroleum Survey Group (EPSG), which is engaged in standardization in the field of coordinate systems and is the holder of the register of their identifiers - EPSG codes - refused to assign its own official code to the Web Mercator system, motivating it with its notorious geometric imperfection. Therefore, in the network you can find links to this system through unofficial codes: EPSG: 900913, EPSG: 102113 and others. However, in 2008, this organization had to surrender and assign a code, since the popularity of the system grew, and it had to be somehow uniquely designated so as not to generate even greater anarchy. The first attempt to define the system was not entirely successful, but in the end it was assigned the official SRID EPSG: 3857.

Algebra

Since projections are the subject of mathematics, I will start with formulas, and then give them a graphic illustration. Strictly speaking, it is not even necessary to be proficient in trigonometry in order to understand the difference between the implementation of coordinate systems based on the Mercator projection, the sphere in one case and the ellipsoid in the other. Formulas vary markedly in appearance.
The projection of the Mercator of the ellipsoid on the plane is given as follows:

x = a × λ
y = a × ln [tan (π / 4 + φ / 2) × ((1-e × sinφ) / (1 + e × sinφ)) ^ (e / 2)]

Where:
x and y - rectangular coordinates,
λ is the longitude on the ellipsoid in radians,
φ - latitude on the ellipsoid in radians,
a - the value of the major semiaxis of the ellipsoid,
e - the value of the eccentricity of the ellipsoid (the relationship of the major and minor semiaxes).

If instead of an ellipsoid a sphere is used, as it happens in the Web Mercator coordinate system, everything becomes much simpler, since the formula for ordinates (Y-axis) degenerates, giving the following:

x = a × λ
y = a × ln [tan (π / 4 + φ / 2)]

Agree, it looks much simpler and shorter, which was what Google developers wanted. This allows quite a noticeable reduction in the number of mathematical operations when working with cartographic materials in client and server applications.

Geometry and Cartography

Even if you do not go into the formulas at all, simple illustrations demonstrate the essence of the problem quite well. Let me explain first that the principle of construction of the Mercator projection is that any point on the surface of an ellipsoid or spheroid is projected onto a cylinder, inside which this ellipsoid is placed so that their vertical axes coincide, and the surfaces either touch along one line (the most frequent case), or intersected by two. (See illustration above). Further, the conditional projection rays exit the center of the ellipsoid, intersect its surface at point P and fall on the surface of the cylinder at point P ' , to which the corresponding point of the Earth's surface is transferred. It is easy to imagine that if the real surface of the Earth is first projected not onto an ellipsoid close enough to its shape, but onto an idealized sphere, then when the points of the sphere are projected onto the cylinder, the same starting points of the earth's surface will be at a different distance from the line equator along the vertical axis than in the case of the ellipsoid.

Let me try to illustrate the scale of the disaster. Take NASA EOSDIS archive satellite image in natural colors of the Central Federal District of Russia, made by MODIS Aqua with a resolution of 250 meters per pixel on September 21, 2014 (this very day - because it was clear, so it will be more beautiful) - this will be our background.

Next, using the Overpass Turbo request, we unload the administrative boundaries of the Moscow region in the GeoJSON format from the OpenStreetMap database. Request Code:

<osm-script output="json" timeout="25"> <union> <query type="relation"> <has-kv k="name" v=" "/> <has-kv k="boundary" v="administrative"/> <bbox-query {{bbox}}/> </query> </union> <print mode="body"/> <recurse type="down"/> <print mode="skeleton" order="quadtile"/> </osm-script> 

Now, using the Global Mapper, we transform the data of the boundaries of the Moscow region from a geographical projection into a Mercator projection of the WGS84 ellipsoid. And further, in order to simulate the situation when the coordinate system is identified incorrectly, copy the resulting data and manually change the definition of the coordinate system on the Web Mercator. In reality, the opposite is possible: the data in Web Mercator can be taken as data in WGS84 / Mercator (this is more than possible, because Web Mercator has a lot of names, some of which have “WGS84”), but our imitation, it will differ only in the direction of the shift. The resulting data will be loaded into the Global Mapper, we will impose a grid over 100 kilometers and see what happened.



The green outline on the map is where it should be, and the red outline is shifted. The magnitude of this shift is 19.6 kilometers. This does not mean that such an error exists in all map services using this coordinate system, not at all. But it will manifest itself in the case, if we take the data in this system and try to combine it with other data without its correct accounting. In this case, an incorrect inverse transformation to geographic coordinates will be applied to it, which will lead to an error.

Navigation

Some cartographic projections have special properties that are critical for solving navigation problems. The Mercator projection is one of them, because it is widely used to create sea and aeronautical charts. This is possible due to the geometrical property of this projection, as conformity. In this case, it means that the shape of objects of sufficiently large size on this map is preserved, since the values ​​of the angles between the lines are preserved. For navigation, this means that looking at the map, you can calculate the direction to the desired point relative to the meridian (the direction to the geographical north) and, moving in this direction on a magnetic compass or at a constant angle to the line to the North Star, be in the right place. Such a path is called a “ loxodrome ” and is not the shortest path between two points on the surface of the Earth, and modern navigation devices allow us to calculate the path from the “ orthodrom ” - really the shortest line, but the Mercator projection is not refused, because the map made in it gives the ability in an emergency situation to use the tools at hand to navigate without relying on the GPS receiver and other electronics.

And here the Web Mercator coordinate system is deceptive. Although it is based on the Mercator projection, the use of a sphere with a constant radius, as the ultimate simplification of the model of the Earth’s surface, deprives its properties of conformity. This means that moving with a constant course angle, measured on such a map, it will not be possible to get to the desired point because of the distortion of angles in this coordinate system. It would seem that this is not so important for web services, because no one in his mind will pave the way for them in an emergency. However, the variety of web services is large, and it is impossible to guarantee that one of the developers would not even consider the directions in this projection. When calculating in this projection, the error can accumulate very much. Plus, nowadays, tools like SAS.Planet are very popular, pumping data from web services, and no one can predict what the user will do next.

The scale of the problem in this case is also quite easy to measure. Take the same picture for the background, the same data on the position of the administrative boundary of the Moscow region. Now we need three lines: the orthodrome (the shortest, taking into account the curvature of the Earth) and two loxodromes, built in the Mercator / WGS84 and Web Mercator systems. We will build these lines between the most southerly point in Serebrynoprudsky district of the region, not far from the settlement with the amusing name “Mochily” and the most northerly - in the Taldomsky area.

Construct the orthodrome. Now we measure its length (a little less than three hundred and seven kilometers) and the initial angle relative to the meridian. Next - the most interesting. We will reproject the workspace into the Mercator projection and construct a line in this projection from the same initial starting point, specifying the measured angle and the length of 307 kilometers, without looking where it will reach with the other end. Repeat the same, but in the coordinate system Web Mercator. Two loxodromes are ready. For clarity, we will also find a central point at the orthodrome, dividing it in half and placing a marker in this place. We re-project the workspace in UTM 37N WGS84 in order to achieve minimal distortion of angles, proportions and other map properties.



On this scale, almost nothing can be disassembled - all the lines almost merge. But let's take a closer look at the center of the lines, including a pre-grid with a step of 100 meters.



The green line with a black dot on the map is the orthodrom. Yellow - Loxodrome, which is built in the Mercator / WGS84, red - Loxodrome in the Web Mercator. As expected, the loxodromes went from the orthodroms, because they are not the shortest distances and relatively straightforward orthodroms are arcs. Thoroughly gone - more than 500 meters. But where did they lead us?



The yellow loxodrome, built in the Mercator projection of the WGS84 ellipsoid, described the correct arc, “magically” returned to the desired point. This means that in this projection one can get to the desired point, knowing the initial course angle and moving all the time at this angle to the direction to the geographical north. And it didn’t happen to the red one - it missed by more than one and a half hundred meters. One and a half hundred to three hundred thousand meters of the way. Four hundredths of a percent. Is it a lot or a little? It is enough not to consider it conformal and not to use for calculations where it is important.

Names, appearance, passwords

The problem with determining whether a Web Mercator coordinate system is being used is not made up. Because of her, let's say, “anarchic” past, she has so many names that everything is simply impossible to enumerate. However, I will try to demonstrate how terrible everything is by listing only some of the known names and codes of this coordinate system:

Web Mercator, Google Web Mercator, Spherical Mercator, WGS 84 Web Mercator, WGS 84 / Pseudo-Mercator (despite the fact that the “pseudo” is not Mercator, but WGS84) , WGS84 Web Mercator (Auxiliary Sphere), Popular Visualization CRS / Mercator, WGS84 / Simple Mercator.

EPSG: 900913, EPSG: 3785, EPSG: 3857, EPSG: 102113, ESPG: 102100, EPSG: 41001.

This is how the system looks in PROJ.4 format:
+proj=merc +lon_0=0 +k=1 +x_0=0 +y_0=0 +a=6378137 +b=6378137 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
Here it is necessary to pay attention to equal values ​​of the parameters of the dimensions of the semi-axes of the ellipsoid a and b . Their equality means using a sphere. In case this is an “honest” Mercator projection of the WGS84 ellipsoid, it is EPSG: 3395 , in the PROJ.4 format it is defined like this:
+proj=merc +lon_0=0 +k=1 +x_0=0 +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m +no_defs

Another unpleasant situation is the existence of definitions of this coordinate system in Well-Known Text (WKT) format, in which the definition itself in WKT almost coincides with the definition for an “honest” coordinate system using the WGS84 ellipsoid, that is, there is the string SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"] , but then either the PARAMETER["semi_minor",6378137.0] declaration overriding this line or EXTENSION["PROJ4","+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +wktext +no_defs"],AUTHORITY["EPSG","3857"] . The problem with such definitions (going as a PRJ file with some data) is that no one guarantees you that the software with which you will open this data will understand what exactly they want from him in such a contradictory definition as the first option. And the introduction of the PROJ.4 format into the definition, which contradicts everything that is written in the WKT, is generally a subtle perversion, since it is not known who generally supports this. It turns out that, at best, the program will generate an error about the unreadable content of the coordinate system definition, and at worst, it will try to interpret this vinaigrette, ignoring the incomprehensible, which will cause the data in the Web Mercator to be read as the data in the “honest” WGS84 / Mercator.

The logic of the situation

I’m not trying to prove here that the Web Mercator coordinate system is no good. It suits, of course. And exactly the same conclusion (except for issues where compliance with US military standards is important) can be found in the NGA report. It is just important to understand the difference between coordinate systems and their capabilities. It is important to understand that Web Mercator is used almost everywhere: Google, OpenStreetMap, Bing, Yahoo and countless other services. It is also embedded in the Slippy Map Tiles format, which stores many tile sources of raster data. It is so popular that not everyone who uses it thinks about how exactly it works. And sometimes it is worth thinking, especially if the planned service should perform functions that are more complex than simply displaying a picture with a map.

Some amusing facts instead of a conclusion

Agency NGA, with a report which began a new round of this story, before the advent of services such as NASA World Wind, Google Maps, Yandex.Maps and others, was the only available source of satellite images of relatively high resolution (10 meters per pixel, black and white image ) to the territory of Russia, which could be downloaded free of charge through the NIMA Raster Roam service (then NGA was also called NIMA - National Imagery and Mapping Agency). These images were part of a reconnaissance program carried out by satellites since the fifties, and were included in the declassification program in 1995.

The Yandex.Maps service does not use the Web Mercator coordinate system, it uses an honest Mercator projection of the WGS84 ellipsoid, EPSG code: 3395. With what it is connected initially, I do not know, but it would be very interesting to hear the comments of Yandex employees, who are present here, at Habré, in considerable quantity.

The local cartographic services of the Scandinavian countries often do not use the Mercator projection at all, preferring the coordinate systems adopted in these countries, for example, the Norwegian public service Norge i Bilder uses three UTM projection zones and EUREF89 datum. This is due to the fact that in northern latitudes the Mercator projection gives too strong scale deformations.