Solved

Mercator to Latitude Longitude Conversion Script

Posted on 2007-11-19
25
2,292 Views
Last Modified: 2012-08-13
Hi,
 This site

http://www.jstott.me.uk/jscoord/

has a great script to convert utm coordinates to lat/lng. However, I need to modify function below (taken from the script which can be download from above), to do mercator - lat/lng conversion.
Thanks for the help!!!
function UTMRefToLatLng() {

  var wgs84 = new RefEll(6378137, 6356752.314);

  var UTM_F0   = 0.9996;

  var a = wgs84.maj;

  var eSquared = wgs84.ecc;

  var ePrimeSquared = eSquared / (1.0 - eSquared);

  var e1 = (1 - Math.sqrt(1 - eSquared)) / (1 + Math.sqrt(1 - eSquared));

  var x = this.easting - 500000.0;;

  var y = this.northing;

  var zoneNumber = this.lngZone;

  var zoneLetter = this.latZone;
 

  var longitudeOrigin = (zoneNumber - 1.0) * 6.0 - 180.0 + 3.0;
 

  // Correct y for southern hemisphere

  if ((ord(zoneLetter) - ord("N")) < 0) {

    y -= 10000000.0;

  }
 

  var m = y / UTM_F0;

  var mu =

    m

      / (a

        * (1.0

          - eSquared / 4.0

          - 3.0 * eSquared * eSquared / 64.0

          - 5.0

            * Math.pow(eSquared, 3.0)

            / 256.0));
 

  var phi1Rad =

    mu

      + (3.0 * e1 / 2.0 - 27.0 * Math.pow(e1, 3.0) / 32.0) * Math.sin(2.0 * mu)

      + (21.0 * e1 * e1 / 16.0 - 55.0 * Math.pow(e1, 4.0) / 32.0)

        * Math.sin(4.0 * mu)

      + (151.0 * Math.pow(e1, 3.0) / 96.0) * Math.sin(6.0 * mu);
 

  var n =

    a

      / Math.sqrt(1.0 - eSquared * Math.sin(phi1Rad) * Math.sin(phi1Rad));

  var t = Math.tan(phi1Rad) * Math.tan(phi1Rad);

  var c = ePrimeSquared * Math.cos(phi1Rad) * Math.cos(phi1Rad);

  var r =

    a

      * (1.0 - eSquared)

      / Math.pow(

        1.0 - eSquared * Math.sin(phi1Rad) * Math.sin(phi1Rad),

        1.5);

  var d = x / (n * UTM_F0);
 

  var latitude = (

    phi1Rad

      - (n * Math.tan(phi1Rad) / r)

        * (d * d / 2.0

          - (5.0

            + (3.0 * t)

            + (10.0 * c)

            - (4.0 * c * c)

            - (9.0 * ePrimeSquared))

            * Math.pow(d, 4.0)

            / 24.0

          + (61.0

            + (90.0 * t)

            + (298.0 * c)

            + (45.0 * t * t)

            - (252.0 * ePrimeSquared)

            - (3.0 * c * c))

            * Math.pow(d, 6.0)

            / 720.0)) * (180.0 / Math.PI);
 

  var longitude = longitudeOrigin + (

    (d

      - (1.0 + 2.0 * t + c) * Math.pow(d, 3.0) / 6.0

      + (5.0

        - (2.0 * c)

        + (28.0 * t)

        - (3.0 * c * c)

        + (8.0 * ePrimeSquared)

        + (24.0 * t * t))

        * Math.pow(d, 5.0)

        / 120.0)

      / Math.cos(phi1Rad)) * (180.0 / Math.PI);
 

  return new LatLng(latitude, longitude);

Open in new window

0
Comment
Question by:chrismarx
  • 14
  • 11
25 Comments
 
LVL 53

Expert Comment

by:Infinity08
ID: 20313107
Doesn't the M in UTM stand for Mercator ?
0
 
LVL 53

Expert Comment

by:Infinity08
ID: 20313119
0
 

Author Comment

by:chrismarx
ID: 20313342
mercator projection is not the same

http://en.wikipedia.org/wiki/Mercator_projection

however, a conversion from UTM - Mercator would also work if anyone has that!!
0
 
LVL 53

Expert Comment

by:Infinity08
ID: 20314085
UTM uses Mercator projections, so you can see where I misinterpreted your question ;)

To convert from the standard Mercator projection (x, y) coordinates to a latitude, longitude pair, do this :

        lat = tan^-1(sinh(y))
        lon = x + lon0

with lon0 the longitude of the center of the map. Just put that in JavaScript ...
0
 

Author Comment

by:chrismarx
ID: 20314571
ooooh, that looks nice, let me just write out what i think the calculation should be

lat = tan raised to the power of (-1 * (sin (h??) (y)))

what is h? is that right?
0
 
LVL 53

Expert Comment

by:Infinity08
ID: 20314809
tan^-1 is the inverse tangent - I think that's Math.atan in JavaScript.
sinh is the hyperbolic sine - I'm not sure if that exists in JavaScript, but if not, you can use this equality :

        sinh(x) = (exp(x) - exp(-x)) / 2
0
 

Author Comment

by:chrismarx
ID: 20314860
thanks, i will give this a try today and get back to you
0
 

Author Comment

by:chrismarx
ID: 20327275
ok,
 below is the javascript, but it does not seem to give the correct results. Im particularly curious about the forumla for lon -- how can lon be the mercator X coords + center lon? a typical mercator coordinate is
(x,y)
1851604, 331835
			var sinh = (Math.exp(coordsY) - Math.exp(-coordsY)) / 2;

  			var degrees = centerLatLng.substring(0,centerLatLng.indexOf("°"))/1;

  			var minutes = centerLatLng.substring(centerLatLng.indexOf("°")+1,centerLatLng.indexOf("."))/ 60;		

   			var seconds = centerLatLng.substring(centerLatLng.indexOf("."),centerLatLng.indexOf("'")) / 3600;

	   		var centerLng = degrees + minutes + seconds;

			var lat = Math.atan(sinh);

        	var lng = coordsX + centerLng;

Open in new window

0
 
LVL 53

Expert Comment

by:Infinity08
ID: 20328444
>> how can lon be the mercator X coords + center lon?

You have to use the same units for the x,y and lat,lon values of course. I would prefer radians. (if they're not the same, you'll have to convert them).

That way, the center longitude of the Mercator map added to the Mercator x coordinate would then obviously give you the longitude.


>> 1851604, 331835

What's the unit ? That doesn't look like standard Mercator coordinates ...


Finally, also check what units the standard JavaScript functions expect ... Probably radians, but I'm not sure.
0
 

Author Comment

by:chrismarx
ID: 20328529
ok, i'm a little rusty here, those coords are in meters. i have no idea how i would convert that to radians, do lat lng also convert to radians?

i think you're right about the radians things, Math.exp(coordsY) returns infinity....
0
 
LVL 53

Expert Comment

by:Infinity08
ID: 20328635
>> Math.exp(coordsY) returns infinity....

That's for sure for such a high value of y lol.


>> ok, i'm a little rusty here, those coords are in meters.

That's not the standard Mercator projection afaik.

Luckily it's easy to convert meters (m) to radians (rad) :

        rad = m / R

with R the radius of the Earth (average here) :

        R = 6372795.477598 m
0
 

Author Comment

by:chrismarx
ID: 20329022
ok! great, see the script below. Lng is now returning the correct coordinates, but lat is still way off, perhaps you can spot the problem. sinh is returning numbers like 0.0345079, the Math.atan() doesn't seem to change this number either-
var rad = 6372795.477598;

var sinh = (Math.exp(coordsY/rad) - Math.exp(-coordsY/rad)) / 2;

alert(sinh);

var lat = Math.atan(sinh);

var degrees = centerLatLng.substring(0,centerLatLng.indexOf("°"))/1;

var minutes = centerLatLng.substring(centerLatLng.indexOf("°")+1,centerLatLng.indexOf("."))/ 60;		

var seconds = centerLatLng.substring(centerLatLng.indexOf("."),centerLatLng.indexOf("'")) / 3600;

var centerLng = degrees + minutes + seconds;

var lng = (coordsX/rad) + (centerLng);

Open in new window

0
Do You Know the 4 Main Threat Actor Types?

Do you know the main threat actor types? Most attackers fall into one of four categories, each with their own favored tactics, techniques, and procedures.

 
LVL 53

Expert Comment

by:Infinity08
ID: 20329773
>> sinh is returning numbers like 0.0345079

Depending on the input, that can be right ... What was the input ? Remember that the result is in radians, not degrees !!


>> the Math.atan() doesn't seem to change this number either-

For low values, that's normal.


Also :

>> var centerLng = degrees + minutes + seconds;

That has to be in radians too !! Convert degrees, minutes and seconds to radians.


For the example input you posted earlier :

>> 1851604, 331835

The result should be (with center longitude = 0) :

        lon = 0.29054816 radians = 16.6471833 degrees
        lat = 0.0520470412 radians = 2.9820758 degrees
0
 

Author Comment

by:chrismarx
ID: 20331216
hmm.... yeah, something is still not working here, i really appreciate your patience. here's all the info:

here are the coords of the center of my map in mercator (this is western ny on my map)
1143031.34, 835295.19

the center longitude is 81.16691694444445 (and the correct center lat should be 42. etc)

the new code is below. none of these values correspond to what im expecting. when i reproject the map im using to utm, and use the above script, i get the correct lat/lng, so i think the mercator coords are correct.




var rad = 6372795.477598;

var sinh = (Math.exp(coordsY/rad) - Math.exp(-coordsY/rad)) / 2;

var lat = Math.atan(sinh)*(180/Math.PI);
 

var lngrad = (coordsX/rad) + ((centerLng)*(Math.PI/180))

var lng = (lngrad)*(180/Math.PI);

Open in new window

0
 
LVL 53

Expert Comment

by:Infinity08
ID: 20331674
>> here are the coords of the center of my map in mercator

??? The center of the map should be the (0,0) point, and it should have a corresponding (lat,lon) pair, usually, the equator and the 0-meridian.


>> the center longitude is 81.16691694444445 (and the correct center lat should be 42. etc)

Are you sure this is a Mercator map ?


>> when i reproject the map im using to utm

Which is it now ... is it UTM as I suspected earlier ? Or is it a basic Mercator projection ?
0
 

Author Comment

by:chrismarx
ID: 20331844
oh, i thought you meant the center of my current map, which, being in western ny, are those coordinates.

yes, i am working with a gis internet mapping application, and it returns enough information for me to calculate the correct projected coodinates; it also returns the center lat/lng (but nothing more).  In the GIS, the projection is indeed mercator. When i reproject the map to use UTM, the coords are very different. are there different versions of mercator? is there any way to check whether the coords im using represent what you think of as mercator?
0
 
LVL 53

Expert Comment

by:Infinity08
ID: 20333326
>> are there different versions of mercator?

A Mercator projection is just a type of projection. The original one is sometimes still used for maps


>> is there any way to check whether the coords im using represent what you think of as mercator?

Is this a world map ? If not, what does the map contain ? What are the maximum and minimum (Mercator) coordinates ? What's the center latitude/longitude ?

Can you also give a known correct pair of Mercator coordinates and the corresponding latitude/longitude ?
0
 

Author Comment

by:chrismarx
ID: 20334858
ok, below is the projection/datum info. i have a series of maps, but they're all generally around ohio/western ny.

mercator
World Geodetic 1984 (WGS84) Auto
Center Latitude 0.00
Center Longitude 0.00

Local offset -9367812.99, 4954053.105
local scale 0.298582119, .29858213
scale correction 1.00, 1.00, 1.00
false easting/northing 0.00, 0.00


here are some coordinate pairs
1059121.38, 855806.52  || -81.3117, 42.3236
1661408.43, 1124318.39 || -79.6963, 42.8539
0
 
LVL 53

Accepted Solution

by:
Infinity08 earned 500 total points
ID: 20335135
>> Local offset -9367812.99, 4954053.105

This is the offset to the center of the map, and with the earlier formulae, you get this corresponding latitude, longitude pair :

        lat0 = atan(sinh(4954053.105/R)) = 0.70916 rad = 40.63195°
        lon0 = -9367812.99/R + 0.0 = -1.46997 rad = -84.22303°

So, that's our center.

Now the coordinate pair examples you mentioned :

>> 1059121.38, 855806.52  || -81.3117, 42.3236

        lat = atan(sinh(855806.52/R)) + 0.70916 = 0.84305 rad = 48.30321°
        lon = 1059121.38/R - 1.46997 = -1.30378 rad = -74.70085°


>> 1661408.43, 1124318.39 || -79.6963, 42.8539

        lat = atan(sinh(1124318.39/R)) + 0.70916 = 0.88468 rad = 50.68823°
        lon = 1661408.43/R - 1.46997 = -1.20927 rad = -69.28588°


As you see, they're not a perfect match, so it's probably not standard Mercator. This seems to confirm it :

>> World Geodetic 1984 (WGS84) Auto
0
 

Author Comment

by:chrismarx
ID: 20335420
ok, yes, i've dealt with this before. I just need to compensate for the local scale and offset. the code below does that. now (finally - thank you!!), im am getting the correct latitudes given the mercator coordinates. longitude is still off, for this region, im getting around -165. the code is below
var rad = 6372795.477598;

var centerLng = (-9367812.99732577/rad); 

var X = ((coordsX*0.2985821196282)-9367812.99732577);

var Y = coordsY*0.2985821374044+(4954053.105471344);

var sinh = (Math.exp(Y/rad) - Math.exp(-Y/rad)) / 2;

var lat = Math.atan(sinh)*(180/Math.PI);

var lngrad = (X/rad) + centerLng;

var lng = (lngrad)*(180/Math.PI);

Open in new window

0
 

Author Comment

by:chrismarx
ID: 20335550
oh, my bad, i was using the center you specified, not 0.0
now it works, thank you!!!
0
 

Author Comment

by:chrismarx
ID: 20335587
ok, i can open another question for this one if you want, but heres definitely the last question. the forumla is working now, but the lat/lng coords are still off by around .4 degrees. any way to make this more accurate?
0
 

Author Comment

by:chrismarx
ID: 20335607
sorry, .05 degrees
0
 
LVL 53

Expert Comment

by:Infinity08
ID: 20335702
>> sorry, .05 degrees

I'm not sure, but I would suspect that that's due to rounding errors, either in the input data, or in the calculations. Or maybe the verification values are a bit off ?
0
 

Author Comment

by:chrismarx
ID: 20699858
oh i found the problem. i was using your ellipsoid value, but i should have been using the ellipsoid for wgs84

http://en.wikipedia.org/wiki/Reference_ellipsoid

now its accurate to 4 decimal places!
0

Featured Post

What Should I Do With This Threat Intelligence?

Are you wondering if you actually need threat intelligence? The answer is yes. We explain the basics for creating useful threat intelligence.

Join & Write a Comment

This article will give core knowledge of JavaScript and will head in to your first JavaScript program. I am Durvesh Naik and I am here to deal with this series of JavaScript. I will teach you JavaScript in part wise , as its quite boring to read big…
This article seeks to propel the full implementation of geothermal power plants in Mexico as a renewable energy source.
The viewer will learn the basics of jQuery, including how to invoke it on a web page. Reference your jQuery libraries: (CODE) Include your new external js/jQuery file: (CODE) Write your first lines of code to setup your site for jQuery.: (CODE)
The viewer will learn the basics of jQuery including how to code hide show and toggles. Reference your jQuery libraries: (CODE) Include your new external js/jQuery file: (CODE) Write your first lines of code to setup your site for jQuery…

760 members asked questions and received personalized solutions in the past 7 days.

Join the community of 500,000 technology professionals and ask your questions.

Join & Ask a Question

Need Help in Real-Time?

Connect with top rated Experts

20 Experts available now in Live!

Get 1:1 Help Now