Sunday, 22 April 2018

coordinate system - Mercator projection: problem with latitude formula



I'm trying to get the position of a coordinate on this map (source Wikipedia): enter image description here


If I'm correct than the formula for the position of the latitude would be:


    R * ln(tan(pi/4 + latitude / 2))

But I can't get it to work. The map is 707 px wide so


    R = 707 / pi / 2

right? I tried filling in a latitude of 51° like this:


    707 / pi / 2 * ln(tan(pi / 4 + 51 / 2)) = 91


And if you go to 91 px above the equator then you're somewhere in Spain, which is too low. I was looking for the location of Belgium instead.



Answer



The formula and its application are basically correct: the length of the Equator is 2 * pi * R and equating that with 707 pixels tells us the scale of the map at the Equator. However, any formula involving trigonometric functions of angles that involves pi expects the angles to be given in radians. Therefore, the Mercator position of latitude 51 degrees must be computed as


707 / (2 * pi) * ln(tan(pi / 4 Radians + 51 / 2 Degrees))
= 707 / (2 * pi) * ln(tan(pi / 4 + (51/2) * pi / 180))
= 707 / (2 * pi) * ln(tan(1.23046))
= 116.8

pixels above the Equator.


If degrees are treated as if they were radians, then pi/4 + 51/2 radians is the same as 42.0848 degrees (plus four whole circles--8*360--which makes no difference), which indeed places the location in Spain rather than Belgium. In this case it seems like a relatively small error of just 117 - 91 = 26 pixels has occurred, which can be misleading, because at other latitudes the result would be incredibly far off or undefined. For instance, using 52 instead of 51 in the formula gives a negative tangent whose logarithm is undefined. Increasing to 54.977871437821385 places you 3786 pixels below the Equator!



No comments:

Post a Comment

arcpy - Changing output name when exporting data driven pages to JPG?

Is there a way to save the output JPG, changing the output file name to the page name, instead of page number? I mean changing the script fo...